{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Another explanation about PCA\n",
    "\n",
    "\n",
    "<img src = 'pca.jpeg' width=\"width\" height=\"height\"/>\n",
    "<sub>photo credit: Raunak Joshi</sub>\n",
    "\n",
    "\n",
    "In this lab, we are going to view another explanation about Principal Component Analysis(PCA). PCA is a statistical technique invented in 1901 by Karl Pearson that uses orthogonal transformations to map a set of variables into a set of linearly uncorrelated variables called Principal Components. \n",
    "\n",
    "PCA is based on the Singular Value Decomposition(SVD) of the Covariance Matrix of the original dataset. The Eigenvectors of such decomposition are used as a rotation matrix.  The Eigenvectors are arranged in the rotation matrix in decreasing order according to its explained variance. This last term is related to the EigenValues of the SVD.\n",
    "\n",
    "PCA is a potent technique with applications ranging from simple space transformation, dimensionality reduction, and mixture separation from spectral information.\n",
    "\n",
    "Follow this lab to view another explanation for PCA. In this case, we are going to use the concept of rotation matrices applied to correlated random data, just as illustrated in the next picture.\n",
    "\n",
    "<img src=GaussianScatterPCA.svg>\n",
    "\n",
    "Source: https://en.wikipedia.org/wiki/Principal_component_analysis\n",
    "\n",
    "As usual, we must import the libraries that will use in this lab."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np                         # Linear algebra library\n",
    "import matplotlib.pyplot as plt            # library for visualization\n",
    "from sklearn.decomposition import PCA      # PCA library\n",
    "import pandas as pd                        # Data frame library\n",
    "import math                                # Library for math functions\n",
    "import random                              # Library for pseudo random numbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To start, let us consider a pair of random variables x, y. Consider the base case when y = n * x. The x and y variables will be perfectly correlated to each other since y is just a scaling of x."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAFoxJREFUeJzt3X9sVfd5x/HP0+vabK6txAW7FOLZBW8sUmibuTEZ3UqbWiK0ClSq1LTQZlIXFHWpmmVaa2Q2Fi1W0/2BsmrZIki7pkqmZGoiBwWXzk3rTksAxRQGS0kKCW0CYUCSrkGRGgR79sc9ZI651z5f33PvPeee90u68v3xlc+ji3mec77n+8PcXQCA/HlHvQMAANQHBQAAcooCAAA5RQEAgJyiAABATlEAACCnKAAAkFMUAADIKQoAAORUU70DmMn8+fO9p6en3mEAQGbs27fvFXdfEKdtqgtAT0+PJicn6x0GAGSGmf0yblu6gAAgpygAAJBTFAAAyCkKAADkFAUAAHKKAgAAOUUBAICcogAAQE6leiIYAIQY3X9Ctz184G3vdbU1a+/wYJ0iSjeuAAA0hM2jhy5J/pJ06uw5DYyM1yGi9OMKAEDmrd++W08+/1rZz0+dPVfDaLKDAgAg05Zv2aXX37xQ7zAyiQIAILOWDY/pNxe83mFkFvcAAGTS8i27Yif/rrbmKkeTTRQAAJkT0u3TZGIUUBkUAACZsmx4LHbyb28p6OjXP1HliLKLewAAMqNnaGfstu0tBR28Y3UVo8k+rgAAZEJI8u9qayb5x0ABAJBqo/tPBCX/lUs66POPiS4gAKm1efSQHtjzYuz2K5d06MGbr61iRI2FKwAAqTS6/0RQ8t+wopvkH4grAACpM9vSDtP94i5G+swFBQBAqoQu7UDynzu6gACkxsDIeOzkbyL5V4oCACAVBrdOBK3aeYzkXzG6gADUHd0+9ZHIFYCZrTaz58zsqJkNzdDuQ2Z2wcw+ncRxAWTf0k07Yyf/DSu6Sf4JqrgAmFlB0j2Srpd0paTPmtmVZdp9Q9IPKj0mgMawbHhM52Ou5rxySYfuXHdVdQPKmSSuAK6RdNTdX3D3c5IekrS2RLsvS3pE0ukEjgkgwy7O7g1Zzpkx/slL4h7AIkkvTXl9XNLA1AZmtkjSpyR9TNKHEjgmgIwKnd3b19mq8dtXVS+gHEuiAFiJ96aX9bslfc3dL5iVaj7ll5ltlLRRkrq7uxMID0BaDIyMB430ufszH9C6Dy6qYkT5lkQBOC7piimvF0t6eVqbfkkPRcl/vqQ1Znbe3Uen/zJ33yZpmyT19/ez1xvQIEK2b2wysY5/DSRRAJ6W1GdmvZJOSLpR0uemNnD33ovPzew7kh4vlfwBNKaBkfHYyX9ewfTsyJoqRwQpgQLg7ufN7FYVR/cUJH3b3Z8xs1uiz++t9BgAsitkgleTieRfQ4lMBHP3MUlj094rmfjd/U+SOCaA9Fu/fbeOnH4jVlvO/GuPmcAAqiLkhi/bN9YHBQBA4kJ28GKYZ/2wGByARIVu30jyrx8KAIDEhJ75M7u3vigAABKxbHhs9kYRzvzTgQIAoCKh6/qwcXt6cBMYwJyFruNP8k8XCgCAOVm6aWfspZyl4lr+LOecLhQAAMF6h3ZesuLjTNjEJZ24BwAgSE9A8p9XMJJ/ilEAAMQWMsyzq62ZpR1SjgIAIJbQCV57hwerGA2SQAEAMKOLwzzjYoJXdnATGEBZods3drU1M8ErQygAAEpav323nnz+tdjtGeOfPRQAAJfYPHooKPmzd282UQAAvE3omT/DPLOLAgDgLaGze0n+2UYBACCJ2b15xDBQAFq6KX7ybzKSf6OgAAA5F9Lt09XWrKNfJ/k3CgoAkGPLhsdiJ/++zlZm9zYYCgCQU8u37Iq9iQsTvBoTN4GBnBndf0K3PXwgdvu+zlaSf4OiAAA5wuxeTEUBAHJicOuEjpx+I3Z7Zvc2PgoAkAPLhsdi9/fPKxjr+OcEBQBocCFLObe3FHTwjtVVjAZpwiggoIEt3UTyR3mJFAAzW21mz5nZUTMbKvH5ejM7GD2eMrP3J3FcAOUNjIzHHuM/r2Ak/xyquACYWUHSPZKul3SlpM+a2ZXTmh2T9BF3Xy7pbyVtq/S4AMob3DqhU2fPxWrbZKLPP6eSuAdwjaSj7v6CJJnZQ5LWSvrZxQbu/tSU9nskLU7guABKWL5ll15/80KstnT75FsSXUCLJL005fXx6L1yvijp+wkcF8AUF/fujZv8+zpbSf45l8QVgJV4r2TPo5l9VMUC8OGyv8xso6SNktTd3Z1AeEDjCx3jzwQvSMlcARyXdMWU14slvTy9kZktl3SfpLXu/mq5X+bu29y93937FyxYkEB4QGMLTf4bVnST/CEpmQLwtKQ+M+s1s2ZJN0raMbWBmXVLelTS59395wkcE4CKSzuEJv87111VxYiQJRV3Abn7eTO7VdIPJBUkfdvdnzGzW6LP75X015LeLekfzUySzrt7f6XHBvJsYGQ89kgfiaUdcClzD9kErrb6+/t9cnKy3mEAqRMyu7fJxCYuOWJm++KeYDMTGMiYkOTPDl6YCQUAyJDQM3928MJMKABARvTS7YOEUQCADFg2PFZ6ck0J7S0Fkj9ioQAAKTe4dSL2Wv4s7YAQ7AcApFTo3r3M7kUoCgCQQqF79zLGH3NBAQBShr17USsUACBFQvbulUj+qAwFAEiJkDH+fZ2tGr99VfWCQS4wCghIgdCN20n+SAIFAKizkI3bm0wM80Ri6AIC6ohF3VBPXAEAdXBx+8a4+jpbSf5IHFcAQI1tHj2kB/a8GLs9E7xQLVwBADUUmvz7OltJ/qgargCAGpnL3r1s34hqogAANbB0006dD9h87xd30d+P6qMLCKiy5Vt2kfyRShQAoIrWb9+t19+8EKvtvIKR/FFTdAEBVTIwMq5TZ8/Fbv/syJoqRgNcigIAVEFIn79JOsaZP+qALiAgYb1D8ZN/X2cryR91QwEAEnJxdm/c+73zCsaibqgruoCABIRO8DLR54/64woAqNDo/hPBs3vp9kEacAUAVIDtG5FlFABgjkJW82SkD9KILiBgDkI2cZlXMJI/UokCAAQa3DoRe5hne0uBm71IrUQKgJmtNrPnzOyomQ2V+NzM7JvR5wfN7OokjgvUWkiff3tLge0bkWoV3wMws4KkeyQNSjou6Wkz2+HuP5vS7HpJfdFjQNI/RT+r5x8GpFeereohkC/u0r9JUku89maS/qZ68aCB9X5EumlH1Q+TxBXANZKOuvsL7n5O0kOS1k5rs1bSd71oj6TLzGxhAscujeSPKjALewBzduwn0v03VP0wSRSARZJemvL6ePReaJvkkPwBZN2xn1T9EEkUgFLnOtNvkcVpU2xottHMJs1s8syZMxUHBwAoLYkCcFzSFVNeL5b08hzaSJLcfZu797t7/4IFCxIIDwBQShIF4GlJfWbWa2bNkm6UNP3uxQ5JX4hGA62Q9Gt3P5nAsUubv6xqvxoAaqL3I1U/RMUFwN3PS7pV0g8kHZb0r+7+jJndYma3RM3GJL0g6aik7ZK+VOlxZ3TrXooA5sRVHO0T+1HvgNGYajQKyNzT+yfc39/vk5OT9Q4DObF++249+fxrsdo2mXT068zuRfqY2T5374/TlrWAAEnLhsf0mwvxT4ZI/mgELAWB3Fu+ZVfs5N9kYuN2NAwKAHJtYGRcr795IVbbrrZmzvzRUOgCQm6FLOe8ckmHHrz52ipGA9QeBQC5E7p944YV3bpz3VVVjAioDwoAciV0B6+VSzpI/mhYFADkxvItu2L390uc+aPxUQCQC0s37Yy9iYvE3r3IBwoAGl5o8meYJ/KCAoCG1ju0M2i5BpI/8oR5AGhYIcnfRPJH/lAA0JBCk/8xkj9yiAKAhjK6/4R6ApJ/k5H8kV/cA0DDCJ3gxexe5B0FAA1hYGRcp86ei9WWLh+giAKAzAuZ4EXyB/4f9wCQaZtHD5H8gTmiACCzRvefCOrzJ/kDb0cXEDIppM9fYow/UAoFAJkTsrTDvILp2ZE11Q0IyCi6gJApIcm/vaVA8gdmQAFAZizfsit28u9qa9bBO1ZXNyAg4ygAyISQoZ5dbc3aOzxY5YiA7OMeAFItdHYvyR+IjysApFZo8p9XMJI/EIACgNQKSf7c8AXC0QWE1Akd48/2jcDcUACQKiHr+Le3FBjpA1SALiCkxrLhsaDtG0n+QGUoAEiFZcNj+s2FeOmf7RuBZFRUAMysw8zGzexI9PPyEm2uMLMfm9lhM3vGzL5SyTHReHqHdsZO/l1tzSzqBiSk0iuAIUlPuHufpCei19Odl/QX7v77klZI+jMzu7LC46IBhG7fyDBPIFmVFoC1ku6Pnt8vad30Bu5+0t1/Gj0/K+mwJIZs5Nzo/hO67eEDsduzqBuQvEoLQJe7n5SKiV5S50yNzaxH0gcl7Z2hzUYzmzSzyTNnzlQYHtLqzwOSf1dbM8kfqIJZh4Ga2Q8lvafER8MhBzKzd0l6RNJt7v56uXbuvk3SNknq7+8PGRSCDAg982dpB6B6Zi0A7v7xcp+Z2SkzW+juJ81soaTTZdq9U8Xk/6C7PzrnaJFpoUs7bFjRrTvXXVXFiIB8q3Qi2A5JN0m6K/r52PQGZmaSviXpsLtvrfB4yKiQ1TwlhnkCtVDpPYC7JA2a2RFJg9Frmdl7zWwsarNS0uclfczMDkQPOnRzpHdoZ9DG7SR/oDYqugJw91clXVfi/ZclrYme/4eK/6+RQyFLO9DfD9QWawGhanqGdsZuS/IHao+lIFAVSzfFT/4rl3SQ/IE6oAAgURdn94bs3fvgzddWNygAJVEAkJjQMf7tLQXO/IE6ogAgEaHJf8OKbpZzBuqMm8CoGBO8gGyiAKAioRO8Vi7pIPkDKUEBwJyFbOIisXcvkDYUAMxJaPJndi+QPtwERrCBkXGSP9AAKAAIMrh1QqfOnovVdl7BSP5AitEFhNhClnZgBy8g/SgAmFXoME+SP5ANFADMaP323Xry+ddit+/rbNX47auqFxCAxFAAUNbAyHjs/n6JCV5A1lAAUFLoME+SP5A9FABcYnDrBBO8gBygAOBtBrdO6MjpN2K1NUnHGOYJZBYFAG8J6fZZuaSDdfyBjGMiGCSFze7t62wl+QMNgCsABK3o2dXWzDBPoEFQAHIuZHYv3T5AY6ELKMdI/kC+UQBy6OLG7XGxcTvQmOgCypnQHby62prZuB1oUFwB5MjAyHhQ8t+wopvkDzQwrgByYv323UHr+rCOP9D4KAA5wPaNAEqhC6jBkfwBlFNRATCzDjMbN7Mj0c/LZ2hbMLP9ZvZ4JcdEfMu37CL5Ayir0iuAIUlPuHufpCei1+V8RdLhCo+HmAa3TgTd8CX5A/lTaQFYK+n+6Pn9ktaVamRmiyV9QtJ9FR4PMazfvjv2ip5s3A7kV6U3gbvc/aQkuftJM+ss0+5uSV+V1Fbh8TCL3qGditvpwxh/IN9mLQBm9kNJ7ynx0XCcA5jZJyWddvd9ZrYqRvuNkjZKUnd3d5xDQMXZvbc9fCB2+/aWAskfyLlZC4C7f7zcZ2Z2yswWRmf/CyWdLtFspaQbzGyNpHmS2s3sAXffUOZ42yRtk6T+/v74dzBzjNm9AOai0nsAOyTdFD2/SdJj0xu4+yZ3X+zuPZJulPSjcskf4XqHdjK7F8CcVFoA7pI0aGZHJA1Gr2Vm7zWzsUqDw8wGt07E7u9vaXqH7v7MB9i4HcBbKroJ7O6vSrquxPsvS1pT4v0JSROVHBNFm0cPBe3d+9yd11c3IACZw1IQGTQwMh60rg8btwMohQKQMSHDPNtbCjp4x+qqxgMgu1gLKEN6ApJ/X2cryR/AjCgAGRC6g1dfZysbtwOYFQUg5UIneHW1NZP8AcTCPYAU2zx6SA/seTF2e878AYTgCiClSP4Aqo0CkFIhyX/Dim6SP4BgdAGlzODWidgTvCTW8QcwdxSAFFm6aafOByx/R/IHUAm6gFJi2fBY7OTf3lIg+QOoGFcAKRAyu5ebvQCSwhVAnS3fsit28l+5pIPkDyAxFIA6i7uW/8olHXrw5murHA2APKEAZEBfZyvJH0DiKAApR58/gGqhANRZe0uh7GdM8AJQTRSAOjt4x+pLisDFYZ5s3wigmhgGmgKs2w+gHrgCAICcogAAQE5RAAAgpygAAJBTFAAAyCkKAADkFAUAAHKKAgAAOWXuAVtQ1ZiZnZH0y3rHUcZ8Sa/UO4gYshKnlJ1YsxKnlJ1YsxKnlP5Yf8fdF8RpmOoCkGZmNunu/fWOYzZZiVPKTqxZiVPKTqxZiVPKVqyzoQsIAHKKAgAAOUUBmLtt9Q4gpqzEKWUn1qzEKWUn1qzEKWUr1hlxDwAAcoorAADIKQpATGbWYWbjZnYk+nl5mXaXmdn3zOxZMztsZjXdzDdunFHbgpntN7PHaxnjlOPPGquZXWFmP46+y2fM7Cs1jG+1mT1nZkfNbKjE52Zm34w+P2hmV9cqthKxzBbr+ijGg2b2lJm9P41xTmn3ITO7YGafrmV802KYNVYzW2VmB6K/zZ/UOsaKuTuPGA9JfydpKHo+JOkbZdrdL+lPo+fNki5LY5zR57dL+hdJj6f1O5W0UNLV0fM2ST+XdGUNYitIel7S+6J/x/+cflxJayR9X5JJWiFpb52+xzix/qGky6Pn19cj1jhxTmn3I0ljkj6d4u/0Mkk/k9Qdve6sR6yVPLgCiG+tisld0c910xuYWbukP5b0LUly93Pu/j81i7Bo1jglycwWS/qEpPtqFFcps8bq7ifd/afR87OSDktaVIPYrpF01N1fcPdzkh6K4p1qraTvetEeSZeZ2cIaxDbdrLG6+1Pu/qvo5R5Ji2scoxTvO5WkL0t6RNLpWgY3TZxYPyfpUXd/UZLcvZ7xzgkFIL4udz8pFZOSpM4Sbd4n6Yykf466Vu4zs9ZaBql4cUrS3ZK+Kul/axVYCXFjlSSZWY+kD0raW/XIikXmpSmvj+vSwhOnTS2ExvFFFa9cam3WOM1skaRPSbq3hnGVEuc7/V1Jl5vZhJntM7Mv1Cy6hLAn8BRm9kNJ7ynx0XDMX9Ek6WpJX3b3vWb29yp2bfxVQiFKqjxOM/ukpNPuvs/MViUZW4ljVfqdXvw971LxrPA2d389idhmO2SJ96YPmYvTphZix2FmH1WxAHy4qhGVFifOuyV9zd0vmJVqXjNxYm2S9AeSrpP0W5J2m9ked/95tYNLCgVgCnf/eLnPzOyUmS1095PRZX6py73jko67+8Uz1O+pWADSFudKSTeY2RpJ8yS1m9kD7r4hhbHKzN6pYvJ/0N0fTTrGMo5LumLK68WSXp5Dm1qIFYeZLVexy+96d3+1RrFNFSfOfkkPRcl/vqQ1Znbe3UdrE+Jb4v77v+Lub0h6w8z+XdL7VbxPlQl0AcW3Q9JN0fObJD02vYG7/7ekl8zs96K3rlPxJlEtxYlzk7svdvceSTdK+lE1kn8Ms8ZqxUzwLUmH3X1rDWN7WlKfmfWaWbOK39OOaW12SPpCNBpohaRfX+zSqrFZYzWzbkmPSvp8Hc9QZ43T3XvdvSf62/yepC/VIflL8f79H5P0R2bWZGa/LWlAxXtU2VHvu9BZeUh6t6QnJB2JfnZE779X0tiUdh+QNCnpoKRRRSMv0hbnlParVL9RQLPGqmJXhUff54HosaZG8a1R8WzueUnD0Xu3SLolem6S7ok+PySpv45/n7PFep+kX035DifTGOe0tt9RnUYBxY1V0l+qeJL3Xyp2T9Yl1rk+mAkMADlFFxAA5BQFAAByigIAADlFAQCAnKIAAEBOUQAAIKcoAACQUxQAAMip/wMYvBpMOxq2LgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 1  # The amount of the correlation\n",
    "x = np.random.uniform(1,2,1000) # Generate 1000 samples from a uniform random variable\n",
    "y = x.copy() * n # Make y = n * x\n",
    "\n",
    "# PCA works better if the data is centered\n",
    "x = x - np.mean(x) # Center x. Remove its mean\n",
    "y = y - np.mean(y) # Center y. Remove its mean\n",
    "\n",
    "data = pd.DataFrame({'x': x, 'y': y}) # Create a data frame with x and y\n",
    "plt.scatter(data.x, data.y) # Plot the original correlated data in blue\n",
    "\n",
    "pca = PCA(n_components=2) # Instantiate a PCA. Choose to get 2 output variables\n",
    "\n",
    "# Create the transformation model for this data. Internally, it gets the rotation \n",
    "# matrix and the explained variance\n",
    "pcaTr = pca.fit(data)\n",
    "\n",
    "rotatedData = pcaTr.transform(data) # Transform the data base on the rotation matrix of pcaTr\n",
    "# # Create a data frame with the new variables. We call these new variables PC1 and PC2\n",
    "dataPCA = pd.DataFrame(data = rotatedData, columns = ['PC1', 'PC2']) \n",
    "\n",
    "# Plot the transformed data in orange\n",
    "plt.scatter(dataPCA.PC1, dataPCA.PC2)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, what is the direction in which the variables point?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Understanding the transformation model pcaTr\n",
    "\n",
    "As mentioned before, a PCA model is composed of a rotation matrix and its corresponding explained variance. In the next module, we will explain the details of the rotation matrices. \n",
    "\n",
    "* `pcaTr.components_` has the rotation matrix \n",
    "* `pcaTr.explained_variance_` has the explained variance of each principal component"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Eigenvectors or principal component: First row must be in the direction of [1, n]\n",
      "[[-0.70710678 -0.70710678]\n",
      " [-0.70710678  0.70710678]]\n",
      "\n",
      "Eigenvalues or explained variance\n",
      "[1.63233804e-01 1.72929288e-33]\n"
     ]
    }
   ],
   "source": [
    "print('Eigenvectors or principal component: First row must be in the direction of [1, n]')\n",
    "print(pcaTr.components_)\n",
    "\n",
    "print()\n",
    "print('Eigenvalues or explained variance')\n",
    "print(pcaTr.explained_variance_)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$cos(45^o) = 0.7071$\n",
    "\n",
    "The rotation matrix is equal to:\n",
    "\n",
    "$$R = \\begin{bmatrix} cos(45^o) & sin(45^o) \\\\ -sin(45^o) & cos(45^o) \\end{bmatrix}$$ \n",
    "\n",
    "And $45^o$ is the same angle that form the variables y = 1 * x.\n",
    "\n",
    "Then, PCA has identified the angle in which point the original variables.\n",
    "\n",
    "And the explained Variance is around [0.166 0]. Remember that the Variance of a uniform random variable x ~ U(1, 2), as our x and y, is equal to:\n",
    "\n",
    "$$Var(x) = \\frac {(2 - 1)^2}{12} = 0.083333$$\n",
    "    \n",
    "Then the explained variance given by the PCA can be interpret as\n",
    "\n",
    "$$[Var(x) + Var(y)  \\ 0] = [0.0833 + 0.0833 \\  0] = [0.166 \\ 0]$$\n",
    "\n",
    "Which means that all the explained variance of our new system is explained by our first principal component. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Correlated Normal Random Variables.\n",
    "\n",
    "Now, we will use a controlled dataset composed of 2 random variables with different variances and with a specific Covariance among them. The only way I know to get such a dataset is, first, create two independent Normal random variables with the desired variances and then combine them using a rotation matrix. In this way, the new resulting variables will be a linear combination of the original random variables and thus be dependent and correlated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "angle:  45.0\n",
      "rotationMatrix\n",
      "[[ 0.70710678  0.70710678]\n",
      " [-0.70710678  0.70710678]]\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnWGMHOd53//PLYfSnlxwjzFRRyudKLgGmdCUeOWhZssCDZnUtCNLvZCWGddsg+aDYKABQpU9hIoMk7KV6ADCpookX4TYCAopKmlTOVORAcoGGbglQNWk72iZFpnYlkRpZVSMxJVj3lLcu3v74W6Ws7PzvjOz887uzOz/90W6vd2Zd3eP/3nmef/P84hSCoQQQorDUL8XQAghxC4UdkIIKRgUdkIIKRgUdkIIKRgUdkIIKRgUdkIIKRgUdkIIKRgUdkIIKRgUdkIIKRgr+nHSD37wg2rt2rX9ODUhhOSWc+fO/aNSak3Y8/oi7GvXrsXZs2f7cWpCCMktIvJ6lOcxFUMIIQWDwk4IIQWDwk4IIQWDwk4IIQWDwk4IIQWDwk4IIQWjL3ZHQgjJK9MzNRw6cQlv1Ru4vVLG5I51mBir9ntZbVDYCSEkItMzNTzy3MtoNBcAALV6A4889zIAZErcmYohhJCIHDpxqSXqLo3mAg6duNSnFQVDYSeEkIi8VW/EerxfUNgJISQit1fKsR7vFxR2QgiJyOSOdSg7pbbHyk4JkzvW9WlFwSQWdhG5U0ROicgrInJBRP7QxsIIISRrTIxV8cTOjahWyhAA1UoZT+zcmKmNU8COK2YewD6l1A9E5J8BOCci31FK/djCsQkhJFNMjFUzJ+R+EkfsSqmfK6V+sPz//wTgFQDZfteEEFJgrObYRWQtgDEAL9k8LiGEkOhYK1ASkQ8AOAZgr1LqFwG/fwjAQwAwOjpq67SEkIKShwrPrGJF2EXEwZKoP6OUei7oOUqppwA8BQDj4+PKxnkJIenhCmut3kBJBAtKodojgc1LhWdWSSzsIiIAvgbgFaXUV5MviRDSb/zCuqCWYrFeCaypwrOb8w5a9G8jYt8K4D8BeFlEZpcf+2Ol1LctHJsQ0geChNUlicB6MYmtzQpPm9F/Xi4QiYVdKfV/AIiFtRBCMkKYgNbqDWydOtm1wAWJ7d4js3js+Qs4cP8G3F4poxawhm4qPG1F/3lKD7HylBDSQZiACpaETeGmwE3P1CIfX3dHcHWuiUeeexnb1q+xVuFpK/rPSwMwgMJOCAkgqHTeRQD43Q9xBc4kqo3mAp596Q00mgsoyVIyIEmFp63+LnlpAAZQ2AkhAXhL5wG0CazO0uYXuOmZGrZOncTd+1/A1qmTbRF9mKi6m7ULSrUi9W7THbqL1NyN+Vh3GXlpAAZw0AYhuaKXm3e60vmtUydD89+6fPTZ19/FqYtXUKs3AiP/IJJu1rqvO3j8AuqNZutxN+3jfY6JyR3r2t4TkM0GYAAjdkJygyuWSXLbNojS4VCXj37mzOXWRSFOMUvSdMfEWBW33dIZx8ZJIeWlARjAiJ2Q3GDb290t7rlMdw46IQ4S85FhB8MrV+CtegNDy4VQfmykO2zkyPPQAAxgxE5IbsjK5l2UdFAcIb46t5QeObx7E77ymXtT63eepxx5UijshOSELAhT1HRQULrGVOzi9YSnle7Iy5AMGzAVQ0hOiLt5l8ZGa9R0UFC6Ztv6NTh2rhZa0Xp6//ZU0h1RUkhp0euKVQo7ITkhjjClVSUZlg6KImDPnLkc2TJpg363AehHxSqFnZAcEXXzLq2NVl2pvwKw6bEXce3GPJoL+oZhpy5eMbphbKeV0hTVqBeMfmx6M8dOSAFJa6PVVJFabzRbou7itxOazp9GvjutNgBxrKf92PSmsBNSQNLaaPVXpEbBK2C685dEUvGE68Qz6K4jDnEuGP3Y9KawE5JTTCX7cRwgpuMEMTFWxen92yO3dK0MO6Hr+spn7k0lLaETTwESFXbFicL74cZhjp2QHBKWO4660ZokB63Lt/t5b66JsS+9iPpcE7dXyti1uYpTF69Y38wMynlP7liHh4/MduT1FZAoxx2nrXA/3DiiAqq80mZ8fFydPXu25+clpCjo+rVUK2Wc3r/d+FqvAOoqPUsioVG0/6IQFQHwuS2jeHxiY6zXmQhaS9kp4YmdG7H3yGzgawTAq1P3WT9fqoItck4pNR72PEbshOSQbjfkdCPv/CwoFRi5+6PiXZurePrM5VhrV0DrNbbE3ZTzrloc2uHST098FCjshOQQXSqgMuwYJxuZRt758VvygtI2x87VMDLstNoCxOGZM5cxftdqK2JoutAd3r3JWNjVrc89y31juHlKSA4J2pBzSoJfXp83WvDiWuxq9Ubr9bqo+Opcs6vZmG6e2wYm54mpK2NWOmbahhE7ITkkKBVw7f35tn7jQGfUHXXD04ubkjFdFHQ7dWVnCI3movZ17oUjaeQb1m5BF11npWOmbRixE5JTXNvh4d2bAKBD1F28gqxrzrX1w6u1hUeu0HWTk1592y3Ys2XU+BxbEfItK27K2ciwE2kjMysdM21DYSckx3hTCTq8gjwxVsWuzdW21IkC8IPL72HXZr0IvlVvYHLHutgpl7fqDTw+sRF7toxqX5u0EtT9DLwXtuuGuwQvphROXH9/lqCwE5JjwjZDgwphgvq1NJoL+NvzP9eKr5urjmuOdoXz8YmNrTuLIJJEyEnaBuiKh7atX5Pr3DuFnZAcYxJEXS9z3WvqjWagcAvQujiMeKpIw/BfVCbGqtpWBEmsh0nSKbqN1VMXr6TSY6ZXcPOUkByj2ww1FSrF3UBVuLlZq6tndIuOvBWl29avwaETl7D3yCxKy4VQI8MOnCFBc/HmgZyS4Nr787h7/wtd+cF172dIJNIxgzZWH9YUNeUl905hJ6QPxPVO654fZ/iGe4xuGmBtnTqJyR3r8J5mg1ZhKcXjrktXCHV1rgmnJKiUHdQbTQwJ0FxQrfx4N211t61fE1gk5Z7Tf8yoo/1sFzX1Ego7IT0mbn+WKM/3CxWAtkKlsOlFYbjnXLUsyKbnuOvRnctt7Vt2SoHPiWs3PHXxSuhzvGmUKJ993GlVWYO9YgjpMaY+L5M71nWItC7K1qVbgvqYCPRecy8rS4IbC/pnjgw7uN5cNF4gqpUy3lredOyWOH1c7t7/QqRzCeKlrvo9eSkI9oohJKOYeoQHRZM6EdWlVIKi5agiaxJ1AKjPNXF49yZjSscVwiQ9z+OkPKKe6/blC04QQY9nuWVAGHTFENJjTMMmgpwYJQk2Iep6iqe5wefaHk/v3250uJgmLQHAsKOXnrCUh99fvm39GuO5vMfsx9CLfkBhJ6TH6LzTpk6LQdKu67ViGi6RFK/gmgZIuDZC3UWpMb8Ip9T5u0rZXDEa1Nvl2Lkadm2utlkW92wZDewN04+hF/2AOXZC+kBQ/taUSzelGqq+/K+uV7g74CJJiuQ1X947LA9tyn87Q4IP3LqiNYDDFdfHnr/Q6hZZKTs4+MCG1jGT9KGPuuYswxw7IRaxLQa6/K3OiWHKaXczPWnt/hdirzko+A7LQ5vy366X3c19P/LcDzsahtUbTUx+43zrXKbPIAjd95YXIe8WK8IuIl8H8CkAbyulPmrjmIRkhSTj4+LgFeRavdHKuR86cSnUrthoLmDvkVkcOnEpknjp7gJMvdW7ubkPsg16uTrXbJ1P1wWyuaha9seSYeKTn159b1nEVsT+VwD+HMD/tHQ8QjJDL1u7uscLGmixa3MVz770hjYX7z7XFa+zr7/ben5JBJ/92J2tiUUmwa1ovOq6zVIXb3S8quxAZMlFUxl2cL25kMj+6G4Im/Yh/BS1JW8UrAi7Uup7IrLWxrEIyRq9aO0aNoe00VzACz/8uVHUvc/1z/lcUApPn7mMV6/8Eq+900Ct3ghMrbjR85AAnqr/SE6VyW+cb6VXvBcGt9oUCm2tBOLgbgjr7jSCLjpFbckbBbpiCAkhbYuc3+mhE+9uxs/5Of3Td1vCaLpGLKqltEzQxKGgVrYHj18winZzQeEDt65oc6pUytEaijlD0rqoxHG1DIq1MYiebZ6KyEMAHgKA0VFz431CeknYxmja5eVx5pD2kuGVKzDzxY+3fjblrHVtBrzU55rG4wWuwRnCn+68p/V9xBkinfe2AEnombArpZ4C8BSwZHfs1XkJMRFlgy3tifRJ7Idp4k9ZHDx+ITBnve/o+UjHcyNlbzOyIU86qOwM4Van1GZ/DPqMo7pa0v7esgztjmSg0W2wBTlM0hIEndOj33hTFtMzNW1UHmXtTmkpneK/kLZnbwQH7t/Q1edssjUOgpD7sWV3fBbAbwD4oIi8CeCAUuprNo5NSJqYNtLi2uO8kagr1v7ioSCyKOr+HulzN+YTHe+2lSswMVbF1qmTRstmN46VQbY16rCyeaqU+qxS6leVUo5S6g6KOskLYRtpUafm+GeP+nuBB/V0cTciM4m6OVGpVm8k3rh1+7iHOVK6cawkGY1XVOiKIQNNWLMqIJrYmDZAg0QmyhDqNDG5UgTd2xJ1uBfQsAtpFMeK35lj6jI5qDDHTgYaf7VnEFHEJk4kOj1Tw76j57UpmGqljLd/0YCmENMK/o1LL1ElfUiW9ge8F4HSkGDBd1HwOlFMhVHuEOmxL72o7RUTlHbR9ZofBFujDkbsZOBx29A+uXtT153/okairjCFVY+mKeouSYPyVWUHhx68t+VNHxl2OgRFAOzaXG1zGLnDo4GbrQCqlTJ2ba7iyPffaEv7uL1i3FSWrte8/xo1KLZGHYzYCVkmiT0uLBJ1RSarnvVuuDrXbPusrr0/35HCcWehetE5VbZOnWyNzfPi7RWjuzNSQKwN66JDYSfEQ1x7nNdmVxl2IFCY84TbI8NOm4WvaHlfN31l2iuI+p5Nz3N/Z+oWuaBUWz/4QYapGEJi4m7erd3/Ah4+MttqBXB1rtkm6sBSteXeI7Ot8vukeV/T5KGsEvU9m57n/i5ss3vQ3TAu+fsrISQFdD1Qgp7ndbOEpand37u2x6AxboKlqssw9mwZxY+//MnIPVayQJxc9+SOdYFTlby9Yvw5+iCKdlfUDRR2MvAEjVvTec+T5MgbzQU8+9Ib+Jejq9o2+xSA+UWldakAS8Lvttx9L0Jfln4xMuwEjqSLcuGcGKvi0KfvxcjwzQtXZXmD1ptaiTJzddDhaDwy8MQZt2Ya9ZaUkWEH9bmm8fjV5U3KKE23ek3ZKQXOK9WN6jPNNo1CWsfNMhyNR3JPr2ZTxunbbdq8S0qU6s5avRGYrug3JieKzYEX/r8Jd47roDX5CoPCTjJJN/0/ur0Q6MQ66JY+yNboFshUK2VsW78GL/zw51Z6p+sIsgT2G9NnbWvgRdDfxLFztUJH6N1CYSeZJG6UF/dC4LcpOkPtFZS6Tb8wr/v0TA1H/u8bXb7r/DL5zfM4ePwC3mt0ttyNc+E0Mcij7uJCYSeZJG6UF+cfvf8i4EbX3sjbL9Z+Iffn3l3CJgkVleaCauX9/RdVWwMvBnnUXVzoiiGZJO5YM13eO+jxx57vHBgBLIm6v8AljmPmC9MvZ3JTsx94/eRei6LfLROHQR51FxdG7CSTxI3yTMMq3J7i7mtN+W9/lK+7Ezh4/ELr92/VGyg7Qx3FSYOON5K2MfBikEfdxYXCTjJJnL4t0zM1Y1Mtb6R9a4RCIK8g6W7z640mJr95vrWRSVHvxHYkPcij7uJCYSeZJUqU56ZKotBoLkQqLvIK0qqyo02vZNGdEpc9W0bx7EtvWJ/i5I+kbVlXB3XUXVyYYye5xna3RL8gSfYs41Y5dq6Gz37sTqvHLIm05dC/MP1yW08d0z4FsQMjdtIXbEVw3Tgi/IMZ/D70QycuYe+R2cwOmbZJo7mAUxevYGTYseK991d+Ts/U8MyZyx3VtLQppgsjdtJz4jhNwugmj+uKuDscYlXZgQC49v48/vqlyx1zS4tOrd7Afff8atcVre6rgtwuh05c0rZIoE0xPSjspOdEGT4ctdvi5I51HdNzwnB7wBzevQnXm4utoc31RjPxVKGRYSdX3Rddjnz/DSx0uWegsPS+T+/f3hGBm8SbNsX0oLCTnhNWaBInop8Yq+JzW0YDR6Pt2TJqHHWXxjSj+lwTt92Svwxnc0Ehia/n6lwz8PvRibcAre8h6kWcRIfCTnpOWKFJlIjey+MTG3F496aOApjHJzYaC2PSaOY1JJJak7CsE/T9BA3GEACf2zLaaudrKy1HbpK/0ILknrBCk25Kx3U2uKDH3Y3bNCh6Xt6/8ezFe8cVtQMj+7+kA4Wd9BxdoQmw1BtdJxxJc7LTMzUcPH6BZf9dUnZKxtTV7ZVy7A6M7P+SDhR20hf8kXTQ0AQvSUvHw44/yITZOgVoXXwPnbgUmGpyc+ZxI3BbnR9JO8yxk0xg2sjstmlU1OMPOosh6aNXp+5rOV7CcuZxI/Cg47H/S3IYsZNUiFuApPuHL4C2RW4ceGuvRwTQafvKkmDr1Mm27/GJnRu1323cCJz9X9KBwk6s0830o25uyaNePKZnahgagCrSbjF5928sqNb34n6PT+zcqL3YdtOBkf1f7MNUDLFOXLsiEP+WPKpNzn0eRd0OYd+jrd7rJBmM2Il1urUrAtFvyXUXj31Hz+PhI7Ntm33MrdslLK0VNQLv1bDyQYTCTqzTrdMhzi25TlzcyNyN4Cnq9rHhWOkmXUeiw1QMsU4vnA5RxKXRXECp6H13LVL2DCEZ0nxs3lYASegmXUeiYyViF5FPAPgfAEoA/lIpNWXjuCSfpOV08N66ryo7cEoSOuxiQanQwpqis/XDq/GDy++FfgavfPmTrf8P8v17bY1JYWFSuiQWdhEpAfgLAP8ewJsAvi8ix5VSP056bJJfbDsd/EJTbzThDAlGhh3U55pa10t1+aKy7+j5gd1Afe2dRsuiqOtjUxLB9Eyt9Z2ZLs42cuMsTEoXUQn/2EXkXwM4qJTasfzzIwCglHpC95rx8XF19uzZROcl2SLtjbBNj70Y2ApgSPR2PackOPTpe9uaTQ1q5O5Wj25bvwbHztW0n8PIsIMD92/QfndBn6N/uEYUbB1n0BCRc0qp8bDn2cixVwG84fn5zeXHyICQdoe+6Zmatr+LyYPdXFA4+/q7rYvOoIo6cHOg97FzNezaXNXuPVyda2q/u+mZGvYdPW8lN05bZLrYyLEH/YV0/HMTkYcAPAQAo6OjFk5LskLaHfqSbKg9feYynj5zOfEa8sDIsIPrzUXjBcwdhWdqIxD03YXVA3STG2dhUnrYEPY3AXin4d4B4C3/k5RSTwF4ClhKxVg4L8kIUQZnuGmayrADpYD3Gs3IKRtuqIXjlAQH7t8A4GZe3NReV5fj9j7HS9gdD3Pj2cKGsH8fwEdE5G4ANQC/C+A/WjguyQmmjTB/LtU7MLlWb+DhI7M4+/q7eHxiozZPHyZCg44/L+7+d+vUSe33MrljHSa/cR5NTS7L/e7CLhIAm3ZlkcTCrpSaF5E/AHACS3bHryulLiReGckNpv4gYZGewlK65KWfvYOfvH2tJSDegpXJHeuw98hsiu8gv+zZMorHJza2fg6zhLaJsMarXnZK2LZ+TaTN5pIIc+MZxIqPXSn1bQDftnEskj9M1riHIwryP7x9reMxN9frDqxm/u5m7/SSCD77sTvb7nT80Xm90cSQoGUJdV0xYbZH1xoZJup0sWSXxHbHbqDdsdh4o0Z2VbRL2RnC9eZi29SpsMh6aLktb2XYwS+vz2vTL8BSEP/q1H24e/8Lxgup98JCekdUuyN7xRCr+HPqFHW7NJqLAJZSVXuPzBp9/C7u7737GzrcTdCwfY0FpXDsXA3jd61mxJ5B2CuGWEV3C18SgQC4bWWp80Wka8JEPQ7e/HtQvx8/7O2SXRixE6vorImLSuHVqfsAAF+YfhnPnLnMnHmGqPqsp/59E5N1kmQPCjuxSpQeII9PbMT4XavbNlvfvfZ+K81AeodpA9RbQGSyTpLswVQMscrkjnVwSu0+OqckHT7nibEqTu/f3hqU/MTOe+DoesX2GUF7S9ssMzLsoFJ2IAAqy3ZHL05JWr+PU8bPodP5ghE7iURQ8RDQaXEE0OFLbC4o7D0y2/Kiu5Y9/+2/zlfdbxSQi7sJ98LoreoF7LRP5tDpfEG7IwklqBOfUxJAoc06V3ZKuNUZiuS+cHH96bJsySP2oM+8eNDuSFokbakb5HQJGnDRaC7E7qDoHoWibp8ojdg4d7SYUNgLjo3ZknQ+9IY0Jj2ZvjvOHS0u+dgRIl1jY7YknQ/p425kuv3Jbc1qNX13nDtaXCjsBcfGbMkoxSokGWt/pYx9R8+jZrENQ5hrhXNHiwuFveDoIjZTJDc9U8PWqZO4e/8L2Dp1EgDwxM6NGBl2Yp/fUuBZeE7/9N2WmNtqwxC2cdrN3wbJBxT2ghMUbTtDgrkb8y3h9o5B0425A4DhlfG2ZKqVMl594j68NnVfVxcF0j3VSjk0T05venHh5mnB8fuPV5UdXLsx37Ik+jfMTHnXOLfobk/vrVMnQwc1ELtEFWd604sLfewDhq40vFop4/T+7dp2re6Ue1PHP7c9rNv3+9i52kAPkO4lbj1AR9EXKRRRfexMxQwQ0zM1rTC/VW9geqaGIU1SfEgEtXrDWBxaEsHh3Ztwev92nLp4xbqo23SL2MLGapL8I3RbAxzevQlP7t4EAHj4yGxHio0MFkzFFAhTsYnbUVFHZdgxTqF3Hzfd3zUXVasgJg1nhUL2+ru7UbI3zRVUvOWnJIJFpVrf09nX38XThu9Hh9sxk5504oURe0HQbXpOz9QwPVMztsktOyUoBSsRtivog+KsGBK0LqQHH9iAQ5++N9JG8YJSre/p4PELGL9rdezov1K+eR7d3si+o+cZuQ8gFPaCYNr0PHTikjHSfmLnRrzXiN7fxYQr6NvWrwl97nBOOiaaWFTocA/NfPHjeHL3pjbhNVFvNPHfjs6iorkgjAw7HZ0vnSHBwQc2tH7W3SEtKNW6wJPBIf//sggAc7GJKS3iyoWNCFuwZKGbnqnh2LlwIVFZbefYJd6qzYmxKmYPfBx7toxGepeLKnh0Xdkp4cD9G3DowXtbVanVShmHHry3LcVi+v5YTTp4MMdeEMIGXOg2TRWWov3JHetChyK7DAEIamL7bz68NP9y69TJSMcpomPGfxH1DxWJs0MwMuzgwP0bOqYaBRH2/bGadLCgsBeA6Zkarr0/3/G4188c9o8+6ig0AFg17ARGl2d+djV0un3RCYqco0wiCmJ45YrIG5/u8/YdPR+4wTwoex5kCQp7zgnqlQ50RntAvH/0JU2/EtcBEkTWHCu9Jkph0OSOda2BI2HU6o1WgVeU4iH3d/6/B1aTDh7MseecoE1ToDPamxir4iufuVdbQu531QSJtGBpU3RQo789W0bb8tz+n3dtXqrcDWrV4DIxVsWeLaORzidAoMvJxMRYta1LZJzxd6Q4MGLPOXE69Ln/uA8ev4D6sgvm1mVniu4C4UUBOHauhl2bq6lVlboVlDZZMSSYX0x+1PG7VuPxiY2BvwvzkftrDIadIcwZxu0FfQ5RBme456OQDzaM2HNONx363p+/KShX55qtSD0KjeYCTl280hYV2ioGrVbK1kV92BnCT/70t7H1w6sTH0vnLJmeqWHf0fNau2lQjYFJ1E2fAzdBSRQo7Dknboc+nd89Tqm+u9l6ev92vDp1X2S/dhin92+33gWyuagwPVPDa+8kF0RXVL1tjTc99iImvxm8d+G+JsrdkIvbs6fKlrokART2nBM3p2pj41NhqUWBSz3G8GoT0zM1vGfpWC7NBRXamTJqzvv2Srkj+q43msYWArcbNpv9eC/IbKlLksAcewGIk1PV+d2rlTLmPO18w3j6zGW89LN3MHdj0Ur6RIDIbpG4uHntoPft3qeMaCycLq6oxom+va8JOvfIsIPhlSsCXS9sqUuSQGEfALwbd6vKDpyStEWZUf3ufv7h7WvW1pimUdIVxYePzHacRwHGPjrAUk+Wgw8sWUcfjnHx2bW5arQg+u2ofrgJSrqFwl5w/G6NeqMJZ0gwMuygPtcMjAR1EWZU4kb/aeJetCbGqto7grCLym233LSOhvWk9/LMmct4+szllhXy1MUrjL5JT6CwF5yg1EFzUWF45QrMfPHjHc93o8SkFaRX55qpWBfjcsuKm9tIYekWHd4ceVDpvjMkmFcK/m0K98davYFj52r0k5OekWjzVEQeFJELIrIoIqFTPUjv6XYSfRL3hRvR9lvUgaU7FLewp9vCWO9nEbRZfejBe0PfLBtxkV6S1BXzIwA7AXzPwlpICugEekjEWMUY5MrIK66odtOaOMiJ4rV6nt6/HRNj1UgXQnrQSa9IJOxKqVeUUgxDMoDXW+0tZ5/csQ5OqdOjHtanOygyNZXU2/af28bNbQdRKTut91EpOxgZdmKX40e5ENKDTnoFc+wFwFTODkCbJggrUY/qynDPn2XcDcsgd4rreEmC157ozob1fuz0oJNeEirsIvJdAB8K+NWjSqlvRT2RiDwE4CEAGB2NVhAyyJjml/oxTU8CljZLdURND5jWE8fb3Q+8zhig0xsOAGNferG1seq1N8bBeyGM8/0RYhtRFlqtisjfAfjvSqmzUZ4/Pj6uzp6N9NSBJKgVb9kpadMCOgeLm4AxfcP+ocpBxw9bT5Z7sFdDRHV6pobJb57vqB51hqRjSpHu9RRw0itE5JxSKtSowpYCGSQsAvdjagQWltf1DlXW5dzD1pPV3LHbdwVA4P4DsPTegloCNBcV9h09b2zBaxogTkg/SRSxi8jvAPgzAGsA1AHMKqV2hL2OEbsZUwT+6tR9HY+bImoguJp0SJbmbPrximHYetw1Vbr0h6eJ6f27+e9qjGKjslPqKDK69v58q/2xl6DPkBAb9CRiV0r9jVLqDqXULUqpfx5F1Ek4cVvxmhqBBf3uyd2btJ7uoJy7KSJXCB7C3EsqZafDseO+/6C7DW/hUFQazQU8c+ZyW3QeJOoAbY2k/9AVk0F07g2Tq8LkYAn6na5tgNvB0Jv8+sD0AAALHUlEQVQ33rZ+TWg/FR3VhK83sWfLqHbwhYtNkY26/qympsjgQGHvAXE32NLq7BelGdi29Ws6rJPHztW6EmXvxSiNiUvHztUwftdq4+cSp7eLDWhrJFnAiismLoOUY4/rcOnlOpwhwQduXdHWDEwXyeuGW3tx8+3+5mJbp06mJq5h+WzdsG8TI8MOrjcXA/PyQc/Vtd4lxDZRc+yM2FPG5CiJKwDdWuvc0W1+YQ5qBqZrS7ugFMpOySiQCsD15iIO797Utq646ZCw83gJO3ZY4VDQuQ/cv6H1Gm86yn/XEaX1LiH9gMKeMt024fKjqy49+/q7xnaw7utMo9u86BwuJZHWCL0FpVApO/jF9WaHs8Z/0ZqeqWEoQrTvUvXdObjn090xRMlnmwqH1v5KGWd+drV1Dm8Pdb9gj9+1mp51kgso7Cmjy/HG3WDTRf7eTUlvK4GoVaHedUzP1PDL6/OBz3NF1Y3cDz6wQRvde2eDmi4qfrwVon7B1KW04uaz/SLvXd+CUm091P3CzcEXJC+wQCllbM2u1EX4fsn0FzKZ7gz86zh04pKx/YD/HKs0Q6wVlgqCHnv+gvGiMjLsoFKO1nAr7mzXKIRZIVlsRPIKI/aUseVwiePu8Iq56XXetIP/dVHOUTF0dDSt1VRotXXqpPZzsh0xh73fbvdCCOk3FPYeYEOQgrztuo1Ab3ol6HUufrtgnItHksKkoDSUqUOl6bNL0qslyvtlsRHJI0zF5ISgVMTntox2pHkEwLb1azpeN9TZkr0jbdOL4Rq6NFTc/jhA8l4t7KFOigoj9hyhi/y9G6gK7ZH49EwNjz1/IbAvDNAekfqtgbYwDc4OWkeUx911JrGSus957PkLgXcfLDYieYXCnmOmZ2p49qU3jBuoYcU5QS0ETIVK3aAbnO1fR1z3kC0r6fXmYsdjI8MOPeoktzAVkwK6MXW2zxHmT48yAMNtIeBPZ2xbvyZyWqZaKeO1gM1Q71rC6MY9FLdZWhC6z2h45QqKOsktFHbL9KpHdxR/epigVsoOTl28EpjOOHXxCnZtriIgNd+GV3yrCYS2GzujDSupraifkCxBYbdMN5uA3RDFn24SVLfIyCRspy5eCXTdeDdib1lx808oidB2426x4W23EfUTkjWYY7dMryJAXU66JNImbkE5du9MT1P7Xt2avRux9UazbZD1rc5Q63xRZ4d2a3V0f58kZdJNi2RCsg4jdsv0KgLctn5NR5qk7JTwlc/cnNOpG7Ixe+DjreeYouyoa240F3Dw+AU88tzLbe6S9+c7NyWD6NVdThBpVLQS0m8YsVumFxHg9Eyto0e6oLOSFAiPaMMqY6O2vA2aJhTVeqhz3/Sqjzp7wJCiQWG3TFpDMrzoepycungl8Plh+WudsAW9F92cTx1RUlC6zo0lCdu6JYQEQWFPgbQjwDh5/C9Mv9zRAXLyG+fx2PMXQouGgM73ouuyeKszFFjkEyWdo7NsRu0KSQhphzn2HBI1jz89UwucNdpcVLg61+zKjqnLSR+4f0Og733uxnzosXU2Sd3jhBAzjNhzSNQ8/qETlyLNKo3bxdB0R3Lw+IW2VM3VuWaowyVo2DWdKYR0DyP2HBLFyTE9U4u1+WjDjjkxVsVtt3TGCiaHS5yNYEJINBixZ4S4BTqmqNnNg8fBlh0zro8/7kYwISQcRuwZwHYbgrB2A/4WvjbTHnF9/CzpJ8Q+FPYMELVAJ2pzMZMoPrl7E776mU2pFeTEbSvAkn5C7MNUTAaIErXGKbvXtRuoVsptValJMaWPoqaVWNJPiH0o7BkgSi/yOEMlelX9arrQxHHYAOkWdBEyaFDYM0AUIY6Ti+5X9Wu3w59Z0k+IXSjsGSCKEMedMJSl6ldCSG/h5mlGmBir4vT+7Ti8exMA4OEjs20bpDaGStiEm56EZBcKe4Yw2R6z1l42axcaQshNmIrJEGF5a9vplW6mFrlw05OQ7JJI2EXkEID7AdwA8FMA/0UpVbexsEGkl3nrJFOLXLjpSUg2SZqK+Q6Ajyql7gHw9wAeSb6kwSWtvHVQYVM/pxYRQtIlkbArpV5USs0v/3gGwB3JlzS4pJG31uXtdQ3C6GohJP/YzLH/PoAjFo83cKSRt9ZF5rqpRXS1EJJ/QoVdRL4L4EMBv3pUKfWt5ec8CmAewDOG4zwE4CEAGB0d7WqxRca/kXl49yYr+WtdBL6gFMpOiaX8hBSQUGFXSv2W6fci8nsAPgXgN5XSzzJTSj0F4CkAGB8f58wzDzY2MnWY+sZM7lhHVwshBSSpK+YTAP4IwL9TSs3ZWdLgYbM834+pXQFdLYQUk6Q59j8HcAuA78jSRPkzSqnPJ17VgJGmzZF+c0IGj0TCrpT6F7YWMsjE7QMTF0bmhAwWbCmQAbJanh91sAchJFuwpUAGyGK6JM0NXUJIulDYM0LW0iVpbugSQtKFqRgSCPutE5JfKOwkEPZbJyS/UNhJIFnd0CWEhMMcOwkkixu6hJBoUNiJlqxt6BJCosFUDCGEFAwKOyGEFAwKOyGEFAwKOyGEFAwKOyGEFIxcuWL8U4ZovyOEkE5yI+xsSkUIIdHITSrG1JSKEELITXIj7GxKRQgh0ciNsLMpFSGERCM3ws6mVIQQEo3cbJ6yKRUhhEQjN8IOsCkVIYREITepGEIIIdGgsBNCSMGgsBNCSMGgsBNCSMGgsBNCSMGgsBNCSMEQpVTvTypyBcDryz9+EMA/9nwR/YPvt7gM0nsF+H77wV1KqTVhT+qLsLctQOSsUmq8r4voIXy/xWWQ3ivA95tlmIohhJCCQWEnhJCCkQVhf6rfC+gxfL/FZZDeK8D3m1n6nmMnhBBilyxE7IQQQiySCWEXkS+LyA9FZFZEXhSR2/u9pjQRkUMicnH5Pf+NiFT6vaa0EJEHReSCiCyKSC4cBd0gIp8QkUsi8hMR2d/v9aSJiHxdRN4WkR/1ey1pIyJ3isgpEXll+e/4D/u9pihkQtgBHFJK3aOU2gTgbwF8sd8LSpnvAPioUuoeAH8P4JE+rydNfgRgJ4Dv9XshaSEiJQB/AeCTAH4dwGdF5Nf7u6pU+SsAn+j3InrEPIB9SqlfA7AFwH/Nw3ebCWFXSv3C8+NtAAqd+FdKvaiUml/+8QyAO/q5njRRSr2ilCr6xPF/BeAnSqmfKaVuAPhfAP5Dn9eUGkqp7wF4t9/r6AVKqZ8rpX6w/P//BOAVAJkfCpGZQRsi8icA/jOA9wBs6/NyesnvAzjS70WQRFQBvOH5+U0AH+vTWkhKiMhaAGMAXurvSsLpmbCLyHcBfCjgV48qpb6llHoUwKMi8giAPwBwoFdrS4Ow97v8nEexdKv3TC/XZpso77XgSMBjhb7rHDRE5AMAjgHY68swZJKeCbtS6rciPvWvAbyAnAt72PsVkd8D8CkAv6ly7jmN8d0WlTcB3On5+Q4Ab/VpLcQyIuJgSdSfUUo91+/1RCETOXYR+YjnxwcAXOzXWnqBiHwCwB8BeEApNdfv9ZDEfB/AR0TkbhFZCeB3ARzv85qIBUREAHwNwCtKqa/2ez1RyUSBkogcA7AOwCKWuj5+XilV6++q0kNEfgLgFgDvLD90Rin1+T4uKTVE5HcA/BmANQDqAGaVUjv6uyr7iMhvA3gSQAnA15VSf9LnJaWGiDwL4Dew1O3w/wE4oJT6Wl8XlRIi8m8B/G8AL2NJnwDgj5VS3+7fqsLJhLATQgixRyZSMYQQQuxBYSeEkIJBYSeEkIJBYSeEkIJBYSeEkIJBYSeEkIJBYSeEkIJBYSeEkILx/wFlKkzI4ru4fQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.lines as mlines\n",
    "import matplotlib.transforms as mtransforms\n",
    "\n",
    "random.seed(100)\n",
    "\n",
    "std1 = 1     # The desired standard deviation of our first random variable\n",
    "std2 = 0.333 # The desired standard deviation of our second random variable\n",
    "\n",
    "x = np.random.normal(0, std1, 1000) # Get 1000 samples from x ~ N(0, std1)\n",
    "y = np.random.normal(0, std2, 1000)  # Get 1000 samples from y ~ N(0, std2)\n",
    "#y = y + np.random.normal(0,1,1000)*noiseLevel * np.sin(0.78)\n",
    "\n",
    "# PCA works better if the data is centered\n",
    "x = x - np.mean(x) # Center x \n",
    "y = y - np.mean(y) # Center y\n",
    "\n",
    "#Define a pair of dependent variables with a desired amount of covariance\n",
    "n = 1 # Magnitude of covariance. \n",
    "angle = np.arctan(1 / n) # Convert the covariance to and angle\n",
    "print('angle: ',  angle * 180 / math.pi)\n",
    "\n",
    "# Create a rotation matrix using the given angle\n",
    "rotationMatrix = np.array([[np.cos(angle), np.sin(angle)],\n",
    "                 [-np.sin(angle), np.cos(angle)]])\n",
    "\n",
    "\n",
    "print('rotationMatrix')\n",
    "print(rotationMatrix)\n",
    "\n",
    "xy = np.concatenate(([x] , [y]), axis=0).T # Create a matrix with columns x and y\n",
    "\n",
    "# Transform the data using the rotation matrix. It correlates the two variables\n",
    "data = np.dot(xy, rotationMatrix) # Return a nD array\n",
    "\n",
    "# Print the rotated data\n",
    "plt.scatter(data[:,0], data[:,1])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us print the original and the resulting transformed system using the result of the PCA in the same plot alongside with the 2 Principal Component vectors in red and blue"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Eigenvectors or principal component: First row must be in the direction of [1, n]\n",
      "[[-0.70377647 -0.71042148]\n",
      " [-0.71042148  0.70377647]]\n",
      "\n",
      "Eigenvalues or explained variance\n",
      "[0.98540269 0.10248788]\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnXmYVNW5r99V1dVNNShNizeGBoTk5GBiRAjEYwI5ScQ4REUcaGdzokYzOsSgmEHRmIByjtM98RpihmscQidiixqDA2qCRxIhoMYTSaIelUZvQGgCdNFdw7p/7NrVu6r2WLWra+jvfR4f6epda6+q3vXbX/3W931Laa0RBEEQGodItScgCIIghIsIuyAIQoMhwi4IgtBgiLALgiA0GCLsgiAIDYYIuyAIQoMhwi4IgtBgiLALgiA0GCLsgiAIDUZTNU46duxYPWnSpGqcWhAEoW5Zv379Nq31/l7HVUXYJ02axLp166pxakEQhLpFKfWGn+PEihEEQWgwRNgFQRAaDBF2QRCEBkOEXRAEocEQYRcEQWgwRNgFQRAajKqkOwqCUN90b+hh6apNbOlNMK4tzoKjpzBveke1pyVkEWEXBCEQ3Rt6uGrFSySSaQB6ehNcteIlABH3GkGsGEEQArF01aacqJskkmmWrtpUpRkJhUjELghCILb0JvIf0JqozhQ/LlQNidgFQQjEuLb44A9a853Vd3LrQ//O+H2bqzcpIQ8RdkEQArHg6CnEY1EimTRLfvO/OX/dg2zfp53Ljz6o2lMTspRtxSilJgB3AQcAGWCZ1vrWcscVBKE2mTe9A5UcYNQF5zHnpaf56RFnM2bpYuZ9ZHy1pyZkCcNjTwGXa63/qJTaB1ivlHpca/3fIYwtCEKtkUhw4rVfhZeehhtv5PMLFlR7RkIBZQu71vpt4O3sv3cppf4MdAAi7ILQaOzaBSecAL/9LdxxB1x0UbVnJNgQalaMUmoSMB34fZjjCoJQA7z7Lhx7LPzxj3D33XDmmdWekeBAaMKulBoF3A9cqrX+h83vLwQuBJg4cWJYpxWEYc+QVIG+/TYcdRT85S9w//1w4onhji+ESijCrpSKYYj6PVrrFXbHaK2XAcsAZs6cqcM4ryDUC90beli08mV6E0kAxrTGuOaEg8sW4CGpAn3jDTjySNiyBR55xPi3UNOEkRWjgB8Df9Za31T+lAShseje0MOCX75AMjMYz+zoS7LgVy8A5QmwWxWo17i+Iv1Nmwwh37ULHn8cPv7xkucqDB1hROyzgHOAl5RSG7OPfVNr/esQxhaEumfpqk15om6STGtfAmylUIx7HKo9vapAfUX6L7wAn/mM8e+nn4Zp02znIA3Aao8wsmLWACqEuQhCQ+Imsj29CSYvfMSXQNqJsRN51aE2eEb6zz0Hn/0sjBplROoHHeQ4B2kAVntI5akgVBgvkdUMCmT3hh7H4+zE2I54LMqCo6e4HuN0s9nSm4AnnzQi9f32g9/9LifqTnOQBmC1hwi7IFSYBUdPIRbx/lLrJZB+mmy1xWMsPvkQz+jZ6WYz/+0NcNxxMGmSIeqTJvmagzQAqy1E2AWhwsyb3sHS+YfSFo95HrulN0H3hh5mLVnN5IWPMGvJ6lwU7xX5A/SnMr7mZPZ7sTL3v5/h+3dfw45/OgieeQbe+96i5znNwc/chKFDaT30mYczZ87U69atG/LzCoJfKr1AOGvJaluPvC0eoz+VybM74rEop8zo4JEX32ZHX9Jz7I62OM8uPMLzOGsK5hkbf8P3Vv2A5ycczFdOX8S3zzjc9vUWeuzm/Px8SxDKRym1Xms90+s4idgFoQBTvHp6E77976DYRczxWBSlsPWw71n7pi9RB/+2yLzpHYxsaeKCP6xg8ar/5Jn3fYTPzV/EtsgIR0to3vQOFp98CB1tcRTGTUREvfaQjTYEoYBycsP9Yo5T+K3gsuUbbY8P8r3aty2iNZ0P3ckl/3UfD0+ZzWUnXE4yathFbjeHedM7RMhrHInYBaGAoVggdLJ6gnrVdlG/V0YMAFrD17/OJf91H12HHMnFcxfkRB3EM693RNgFoYBKLxC6WT12Fo1TPo0CTpnREdwWSafhC1+AW27h1TPOY9Hcy8hEBs/p++Yg1CxixQhCAQuOnmK7QGgVu3IWV92sHnPR0zr2pw/an7vXvlk0jgaeemWrr4XSHAMDcM450NUF3/kO77/2Wr6/cUvFFoqlSrU6iLALQgFO/rf5eLnVl06WTk9vgllLVufOefNp05g3vYPuDT22wu42li2JBJx6Kvz613DjjXQfeSZLb3iq6HxhIVWq1UOEXRBscFsgLHdx1a3Hi/m4VQTdipZ820MFG2R0H3Z82aLrFY0PxSK0YI947IIQkHIXV+18dDtMEXQb15cXvn270aFxzRpjg4yLLiq7NYCflFCpUq0eIuyCEJByF1fNXPCo8m4zYEbDdrTFY96R7zvvwCc/CRs3GhtkZHc9crOD/ODnxiBVqtVDhF0QXLAr73cqLjKjZ6eWAFbmTe8g46Pqu6015ni+RXMPdn/yG2/AJz4Br71mbJBh2fXISVxVdv5e+InGvd4noXKIsAuCA052A+BYfRmkatVP5LqjL8mlyzcyIhahLR7zn9b4l78Yor51Kzz+ON37fTDvZvPpg/a3TaPUuHv6XnO3Pi5VqtVDesUIggNO/VzserGYC4lOVkZbPMbGa44qek5hWqUbCjjr8IlcP+8Q9wNfeMHYn1RreOwxuvX+tumbTudVwOtLjnM9hfSMqQ7SK0YQysTv4p81SneiN5HMi9rNG4FfUQcjmr577Zt8u/sl54Oeew4+9SlobjYyYKZNc/TDnTx+P98kJBqvbSTdURAccEpL1BjRvJne51egzTS/oJF6IfesfZOZB7YXi+iTTxo++gEHwBNP5HqpO92g0loXRe7xWJRPH7R/Xj69U1GR9IypXSRiFwQH3NISrd653/S9nmyv9aCReiG2PvhDDzlukOEUgZtRtjXqPmVGB/ev76loZ0uh8ojHLggueHnnHVnR9Jsm6OZtB+V/TB/8vvuMNgEf+Qg8+qixpZ2FIH54kHUFYegRj10QQmDe9A7X9LwtvQnfBUdgeNs+0tcBGNnsPGYuLXHZMjjrLJg1y7BfCkQdjNdwyoyOnKceVYpTZtjbKFJU1BiIsAuCC2a068S4tnhOOH3qNX6/JPcNpDn78In2YwBvffu7cNFFcMwxRqS+7762x3Zv6OH+9T2ksydOa83963sCpWBGlHLNyxdqCxF2QXDBzQ+3Fts89crWQJth+GFcW9w+tVFrLvvdPXzt1z+E+fOhuxtaWx3HCdI+wOnbR1pr8dzrCBF2QXDBzYKwetR+rQoFtMa8P3bWm0aHNYrWmu+szu56NOMYw19vbnYdK4i9UpjGaJcSGaSnjFAdRNgFwQW3jBKrR+23/4kGmpvc/XgzO+Xah15m0sJHcouZkUyaJb/535y/7kF+OuMEvjbny8xa+oxn9Ow0NzNts/D586Z38OzCI3h9yXGObQ/Ec69tJI9dGDb42fSh8JhPH7Q/96/v8dx0o28g5XsevQnnTak7sudc/vxbJNODohpLJ7n54Zs4/pXfcevHT+fm2WeBUr7a7Tpt1AGD1sq6N7bz1Ctbbbfqs8uSkUZetY2kOwrDAj8pf07HnDKjI0/0Pn3Q/rmfR8dj7BlI5YlwuSjyN69uSfbzf7oXc8Rr6/j+pz7Psn85peg5bumITimMbuc03xtAWgfUEH7THSViF4YFTguIi1a+nIvQI0rlMkesx1i3nysUf7fou1SsMxjZ38eP77+Ow956mW8e/RXunXas7XPcrBE/tknhbcltqz7Z3q72EWEXhgVO4tabSObEuVDUTazRbrlVo0EYndjF//3lNXz4nb9x6QmXs/JDn3I81s0acduxyQ3zPZPWAfWHCLswLChV3GCwGGje9I5QFw3dqlb3372Dny//NpN39PClk77J4x843HV+hUVU1rWC0fEYsahytIsKbRgT8dHrF8mKEYYFQapDC7H2ZglT7BYcPcV2Xu/fs41f3XslE3e+w3mnLvIU9bMOn5gXURf2hO9NJEHDmFajn3tbPJb7d0dbnLMOnygbYjQYErELwwJT+Kxecd9Aih19/jxys4HXgqOnFC0mxiKK5qYIewaCWTSFGz9v6U3wL8lt/ORX3yTVt5OzO6/nj+M/WPQ8M8LuaIszab849/3+Le5e+yZRpTjjXybw1Ctbi+yiZEbT2tzEhquPKhoPYOaB7eKjNxCSFSPUFX5SFoOMVSjSTrYE5GeK2M1h0sJHfJ9bKXh9ccFmFpYNMv7t9Ot5utX5dbnN0w2zcViY76MwdAxpVoxS6ifA8cDftdYfDmNMQSikUIj95HC7YT7n2odezkXuI2IR+lMZMjaqmUimuaxrIzd3TrNNLewI4OMXxVNr18Kxx8KoUfD448xL7MPvXXq2lyLqZhVp2O+jUHuEZcX8DPhP4K6QxhOEItx6npQjSHuTGct4GZcjDUFe8MsX+OW6N1n72g7SWucsEDubxom8NgGrV5M6/gTeibdx+vHXsnv5myhF6Nk3ZtZPpd5HoXYIRdi11r9VSk0KYyxBcCLMlrJefdbdSGY0z766PfdzWmvuXvsmXc+/xYCPQqW8hcmHHiJ96nxeHX0AZ3d+l62j2qECufEweDOR1ryNj2TFCHWDU0ZK0EwVP3uUloIfUTf7wCxdtYmL515Bat48Xh57IKedsdgQdQeiSjk25Spk1vvbXbNcwnofhdplyIRdKXWhUmqdUmrd1q1bh+q0Qp3SvaGHWUtW5/UAt0sNtE3Le7ELbv4wLGoz/v9iV96vh7LIyEpHdpHy/vU9/OvTD3DLQ//O+o4PcWbn9fTG7Xupm2S05vUlx/EfnYcSi9iLe0TB2YdP5J4vfMx1o2nf76NQt4SWFZO1Yh72s3gqWTFl8mIXPHkd7NwMo8fDnKthame1ZxUabn1dwKO8/cUueOhiSFqi8VgcTrgt9x4FyV4JEwWMjsc49ZnlfPupn/DU+2bwpXlXsTc2wvO5bfEYMNjCwJoVM6Y1xjUnHBzIH5esmPpEesU0KoXCtfMt42doGHF3Wty7vOsF/qPzUPe9N5+8Ll/Uwfj5yety70/U0hNmbmQNVzR1MU5tY4sey42pTlZmZof6ekxGj2ji84/9zOilPmU2l51wOclozPN5EQX/2JvMy9TRQCyqWHrqob4E2U7IZQ/TxiUUK0YpdR/wHDBFKbVZKXV+GOMKNrgJV4PgtIiX1tp1957uDT1kdm62H9TyuFXUl8TuZHxkGxEF4yPbWBK7k7mRNeW9ABuU1nztkf/DJf91H12HHMnFcxf4EnUwRNwu/TKZ1r42vCisRJVdkBqfUIRda32G1vq9WuuY1nq81vrHYYwr2OBDuKqKh7/t5xi3RTyn3XtM8dqSKd7MGaAvfkDuOHMB8oqmLlrVQN5xrWqAK5ps5pxlbmQNa5ov5rWWM1nTfLGvm0Akk2axZYOMK4+9mEzEf3sDN7fUTyZLkK3xhMZAsmLqjdHjgz0+lJg20c63AD1oE1mF28cxXn1d7MTMFK8bU5306fyt4vp0MzcmT8uJvxmxj1PbbMcfp961fbyUCD+WTnLbQ//O6S8+xq0fP51r51yIVsE+dm6ZME43Qevis1P2j6Q3Ni7isdcbc662Xxycc3X15mTiw9/2c4zpGV/e9YJtK107MTNFamVmNiTJ+ubvskXvZ/jm/YcRLRhvix7LeBtx36Lto363CH/lQLEv72eDDC/MjT6W/+EtkhmdtybwNmO5pe90Ji00XntbPMaiuQcDxZtj2CHpjY2LCHu9YRXIWsuK8WMT+bSSTHG3y47JpeVZsoOeGzGW7w/MZ2VmtvFfgdAqivut35jqZEnszjyx7tPN3Jiyfy+DRPh+N8iAwcXcwv4vCjhlRgfXzzuEmQe28/sH7+A7enC+HWzjWpbRH8mwMjOb3kSSBb98gVEjmjxFXdIbGxsR9npkamdtCHkho8dnLRabx03iYyCxvfiY+Jiih+w6MubS8gqygw5gK7fGbmcRd7EoeW5eZotdwywz8h3BACkdIUqGNBHiZD32FEXZMX4i/LmRNVyZvI9xXW/B2xnumnsC937QWdRh8IZTOEcNPPXK1tx7Me/pB2Cn+zeGZEZ7dqzskPTGhkeEXQiPCthEjrv32Fg6SkE7u1kSuxOShjCPaY0VCZ3plZuRb4QMWkOTMvrEjFfb8sYw8Yrw50bWsKR/Ga1398K7GeiMM3/Kc6xPHlxyCmWeD+7wbcfuG4NTGqcZqYuoNzayeFppXuyCGybDotHGfzdMts8U8TuWW8aJn4yUMM7jxNROoxBo9ARAGf+3FAYBkNhh/1ynxx3QLllArWqAW2O3s6b5Yj6x9ykK1x7tvPLCY+yyY1ZmZrMweQGbM2PJaMXmzFgWJi/IifZVu++h9We9sCMDZ7bCQTHPLBsv8nxwhwXywjWBM0as5QaHRV7JhhkeSMReSV7sgge/AmmLiCS2Q/eXjX8HsVO8CpPCKlwqY5y/7/k7T6o0Z1z2J+eD/Ng1uFdGdm/o4aN6PzocPG8whHq82sYtsduZkf4L16TOy/3OySsvxC4StvPvASZv7+GA5ZuhX8M5rTBh8KPllGXjRZEPbvONqHBNIBZRfCf+S+IJZ8tGsmEaH4nYK8mT1+WLukkmGbygyKswKazCpTLGWbJmCWeuOJObnrvJ+aA5Vxv2jJUCu+b5lT/ko93/yu8SJ/G75ouZ8Y/H8wpqlq7axA3J4rRGOyIKzok+kZeSuEWP9XyecZx9dkwhH/z7a3TdcyWZlIJ/G5kn6kHGsWL2dwEGe+b8eizPH3Jt9hsRZFSEeFaw50bW0BaPsXT+obQm3rEfU21jTfPFvDrirPK+0Qk1jwh7JXErGgpaUOSVTRJW4VIZ4yw5cgnzPzSfyx+7nGueugbbPkReds2LXXz4j9+hQ+XbCAv1jzj8wU/CojaW930BgIXJC9itW1wLeMAQd6sdYpfrXjiGW3aMlek9r/CLe6+iKZpm17/th35Pfv6933FM4rEot5w2uJFHYcXouc8fyPPv/xrE4kR0BoXxHt3WfDsb9XzmPfihYl/JfI3ZYyNONQZCwyBWTCVxsh3M34UxljmOT4vDlRe7QEVA26TKeYwz2N/8bP7XyF1c99vr+Ef/P7jp6JtQhULjltXz5HXE6c97qFUNcE70iVwUMj6yjVtjt7Ndj6JfxxgV6S8epwCrHWKX6/5kZhpzIhvzc989Fjw/9sYL3Hn/d+kbOYJR5yriYwbnoTVs16O4NnWur4XTuZE1fLP5l7yHbainx0P0apauGmtbMTrhj0sBFztFZ9AY2UAmGW3c4PIorDEQGgbZ87SSPPx1WGfTXSESg3m3l+d9Q37XQjs/PxKFltHGwmRhvnthh8gPHAUv3FtswxSex4bCboyaDGNbvssfI89zvo7xw30/QPTIa+hOz/LuKLiojSAbv2ntGKDmsTkzltkDt/ke14s5f/s9t3cv4fUx43jPOXsZs2/x+7Zdj6JPj/BsMFaYpQNALM4lez6PhqLslltitxeLtA0pHSGC5m32Y5zahv1TFCzq9fuyg9HgXUirgXR3rDYvdhlCWUhsJJxwS/AL3E9hUuFNOpMezBm3LoJC8QLpup9gK6gq6irqUNyL5Lqmn3E2r3CdbuY6NcCuf/yFn3Z/lTWpL9Az8HEAZvzjcT7a/QX0g++isq+lOz2LwxnLAfjv1+9H1IPaIV6c8N/PcNMjN/Hye97HLzo/y+J9fmp73Bh20x7ZDdinUJopiR1qW/HrSCa4tvnnNOv+nOCbY/xD7UMbuzznGSXD5H7jGlw74hL797VSrSiGQRfSWkaEvVLYLUICtLaXfmF7WBhkPLZUsy6CFs3NIUrWGc/5WrMs5kbWGLaJUlzLCPbRigWqn92ZHdzCvQxEMiyK3cUYdg+K2c63SD34NdYkL2B1an5R9FpoK/hFa+jxiJSDtuw9feNv+P6qH/D8hIM5/5Sr+c0+Cx1vLk4plCsHZttH6QWMZpftGLFIFDJ2ZVf5pC1LaIsH5nPryJ+WV2MQJAL3015CqBgi7JUizC6MdrbJXx/L/4D5HTfo+c2IzuVDPa4tnms0dUVTV55N8A1a2Fcrvshe9kTeYnlsGe0qVXSapvReLuUXzM7cVuR/d6ht7squIsYNqIAe7Wy/FAqrU1GSlQv+sKJogwynlEsni6hDbeO1ljPJEMkVRAUllvGXrhhhcPx1+34GPju9dGskaARe611IGxwR9koRxmIm2H+grL69+QFzKtV3Or/Tom4eyjjuhskwsHvQvy/4UC84ekrOY7fLEb+QZvZBcQ4J5qpeHmUk7TZKXZzvbUSkOxhFO7ttZ5jQzfxaHcGxerXvni/g0bI3VeBrJ+fz/t++mbdBxrGx35dUeKSUcY+yiq7jsYFHz8dM68zlw089ovRoOWgEHtb1L5SECHul8Cqv9/u11snSsZJMQFPcGN/tWOv5C+dWhOWrvt0NI5mAB74IwLzpxryXrtrElj77fiqn6xgjFXSS4JPs4XFaOaAg23aL3s82ku7XUQZ0E82WSN/MOnk4czhz9DpGqMGeL272i4lTkVJHNnLPnZ+t3PrUzajfD7B26iFcfPQCjm96ztNGAePm4nVMpUjQwurMNNaOuCQv06ZkYQ8agQdtLyELraEieeyVwi1f20/fchO/X10TO4rPN/N8+/Obc3PFR2aKTufmPW96B88uPILxpy4uKkDSGJHqXGI8Qiuvk+ET9PGGJWo1I2y7SLpFpYmSYrsehdaQ1kYs2652c270iVzpvGlttKld3Bq7nb+1nO24IUYvo2xfUprI4PkzGh7ai/r9ABzWzNQTt3B803O2c7TDbD+gtftmGaEQb8/7W8dnns258Wc5gK0o6zX28NdLa0sRdB8AP+0lrOf0+3kQfCHpjpXCLQK5+cMOX1MnQGE5vtOxdoyeYB/pOM0lyNhuqKjhcZsdGhPbs4+lId6OTmzPsxWeI8Vn6WMfFI/rVsbq0bl879daznRM5UtpiOIvE6aQPt2c6+syN7KGpbEf0qLyc8QHdBNNpIzzpzU8kICXU/CvzfCpFlCK7XoUbez2TDfcrkfxkf5lzI2s4dbY7SXNORCTPwnbXxv8Gw/scbDmChZdC1NmnVJqwXOT8JIJ8nkY5vhNd5SIvVTcGmV5RSBBvtbOuRrfbmvQHYs+cJS/cb3QaWPsxPZBMdHpXOReOPuP0cRTjGQv8K+qj/WWI9zK/ZtUaaIOhn9+U+wOXms5k5tidxSJOsAuPYI+WiCpYXlW1I9sgU+PyJ14DLvZoe2jfROt4drUuTlbqeKiDvD6b/P/xo7rLQWBnN+2FEEi8KDIQmvoSMTuhlOk+2KX0cjLml5oLTryikCcfh9vhytfL3580ehg8463Q/NIY96OlaRGv5FQIvYS2USaI+ljN5oH9GhWJL/EjMhfsumSQz8frYF+jfpFH7yRhuNGwMzifjQZbdysnAQ7o+F9/ffyx5YLaVf2i77m+axjlJrWWT7ZIiXH4rAKFjGBROwBkIi9XNwi3UevLM4ZzyThoUvd7Q0zAnGKlPt32fuKpgj7JbF9cN52om7OpcoR0RSirGEkY4lwvNrJibHbOLepOqIOoPZq1N174M00nBy3FXUwSvPdovAMEV5vOZMxDpk8VlI6QkYb/zcm4b7JdXlxmEvC/Q2TcVxXqXQmi4/GcEIwRNidcPta6vQ1N7nHPQIePd65IhXyuz5arZ6BPcHn78Xo8TWRenYgEX5HK5OJcJzq42E8iqwqxe4M/GwPvGNskMEhsZKGMTfsUD5sI5Vd8FUM/t+4ETs/MbCto6IMLqafVyygYKyPOF3TXpksYfT/r6TNM0yRdEcnwvb9zA+IV/rizreKF7H85KeXMheAFRcSpDdLJTiACM8wkmPYw0kk+DlwBqUJaxC0Bo1C9aZRd/fBruwGGe8L/rEwbZVS/PTi52hftoyvPjk6k2+jTDzcSFN1+iZnxWkxHsJvGVCr2z3WKRKxO+GW3tU8MuBYlgjE88agDKvHK3c9KCpCXjQE2W8H4Yt6KXZBO4onGclsopxFgmUMTf739945zYjU92SMDTJsRN1PumLoC6Ta+MtoyrRfVCQ/kp7aaVula/NEw992EvUHvhhO/3+hIoiwO+Hk+33gKEgFFB2rleJpf+jwI3QwCpgW9Q4uRj34laounNqxD4pHaeVYmrhI7eU/8G7HWw76nTTfvudHqBS2G2SYlBqJl4NZoarAMXRP+/n4WmoNcvix4JyOMSN1t7UboeqIsDvh5Pv99TH7Zluxkfb+JQxuh/diV7D0RTvi7aU9L7ln8MP96JX2OzuFhFKlfw+Io3iAOPN1E99Q/VzD3qwxUR5FUe/mFJG79hgW9Odb4QD3RctqYne1aA0pFAPah21UGEnbBS1W3Hx1LyuxBtZtaoqw1iECIsLuxtROI8I1I103KyXZZwi/U1aDuTA6tZOy7I/EdkPcIyUIkfnhrsQ3ggLKCXCbUdyr43xex7hODfB1+ssW97yI+/UU3NUHcQWfHwlja1fUnVAKRqg0Gu3PqrFet1M74dAzsf0rxdvdFy7dInLJZMmnihW1sngaFLfmRlM7s4uRDuzcHM4ftVRhDuv8WUxBqYRN0aQUd2bb/t6iBtilNT9kBNFyM703JeGXCdgvAme3wj6Vi238bgJSTv56i0r7u+UVRtJ/fQzbAKN55GCtxpPXGde6WUU8eoJzszkfffuHHVVsXSwRe1C8cm7dvoqOHm/YIE6oCv854mMMSygklIIkTbneLWETQXELLVytm/mxSnImCQbKidz/lISuBLwnAp+rrKiD/xteue+e5/PtImm3rK+8SJNBP33nW0aXz0hBxlK0GUaMhhVfgGvbjYI62Sy7qhW1IuxB8cq5nXO1caEXEokZv3OLtitaBayMc3ttxhGQZpWiTzdXbOoKY8OOpbqFLpXiJBIkShH39QNwfwImROHckdA6TC59p5xwt6wvNx89PQAt+wxe//F247q1tpIAaeQFwRunhcgwubpDxs57t/7uxB/kL3LG2/3tcWo20fJNkFg+S9SVAAAgAElEQVSvcjeNUaq/4lkj36CFH+oRPEqKY+ljV5DX81w/PLwX/qkJzmqFliqVtoaA/1cdgZN/NJgFVbiA5/bN0yuiTOwYvP6bRzoHC0HTH6u00FgxqlhRW3+9YsLs21yNHtA3THaO2uPtkEoUd9Abfxi8/kzx8ZM/Ce+8NCSLobXCfSQ5hwQziDhu2JFDa3imH54ZgA81GW0CovUr6oU4evgqCifd4a9j46NXDl4/8XY49oZBb92L0RMCHOfxGfParL1eCVljGrNXTJirzNVasT72BuffJbZbbB6MD2gyAf+zxv747a8ZTcNmnk+12kcNNWcQYwVxXiDDJ9nDO047EWkNj2VFfVoMTmksUYdsWqldXDbj3/I3P3crJEpZfpfYbnwGPnCUezqkSZBduLw+Y17zrFfcvt1XkFCEXSl1jFJqk1Lqb0qphWGMaUuYf/xqXUhTO11y0bPCY36FM/1Kr2KQ42+Ck5eVnuNeZ7ht2AHkNshgrbFBBnNHULXOYhXGNmJf9+NBK8NxAe8t58/Auh8bBW1lX082G247fcakdW+olC3sSqko8APgWOBDwBlKqQ+VO64tYf7xq3khHXsDDmUnxkXvt6WAdRFmaqcRvZ/8o/yF3VjA9gd1whyaeJxWtpHhE+xhE9mbX1rDigRsSBobZBzTMvRlo7WAdS9cW5R7xJ3YbkTzJ/+IQN8GzTqO0RNwXBGw+4xVcaGxEQkjYj8M+JvW+jWt9QDwC+DEEMYtJsw/vuNzdOUXbtyKlFw3SbDgtAgztdN4fPR44wPU1EK9OW5+sW7Y8Qn6eCGVctwgo24opfDMiWTCWOi0RXu2CM5F174/Xwqu2Q6Ldhq2g1O7abvxpHVvqITxie8ArLf+zdnHwifMP75bWfVQ+O1OF73Xh818rt/9IxPbgQyN6sFPI8rvaKUF+FSmj+cSKWODjFkt1Z5aaWTc2/YGxyU5wrLLlSM7N3u3IDApFGw/n1czE2bFhRb7p6BZXSNlygwRYQi7g6dQcJBSFyql1iml1m3durW0M4XZtzlvLBsq7bc7fVi82qmaHwzzNRemiDnaONVtzVtJphBlze447+mDns80O26QUT8M0d9q9ATn1gImppXT5CHsdgGW1+fVLghJJYz1IjNNUza5Lomy0x2VUh8DFmmtj87+fBWA1nqx03Nqbms8ry3BKpUWae7GFDRd0dwyzC5FbBgzkM7QHG1M26kizDzfaC3g5rWrqGEPWZvGxeLGDeGvj+V/JsC+DYHT56XULSSH8ZZ5Q5nu+DzwAaXUZKVUM3A6sDKEcYcON+/eKy2ynKKKqZ3Be7vDYNm3XU/sYUwpou6n13rD8sK93imLOl3cCTSZMETdmsYHzm0IVlxo32bAK4FBMmVKpmxh11qngK8Cq4A/A11a65fLHbdsggiumxfolhYZRi58KRdprNW9J7YQCLv11YYS+2gztnZLOUFB4XXr2s43+2YWfj68kiG8Ai7x3h0J5Xur1vrXWut/1lq/X2v9vTDGLIugguvmBbrlAdv52UG9edeMAwfvs7A6VSgZp6SZekumKcJ6LZ/4gwqMX3Dd+t20xfr58FpcddvsRrx3VxrTkCyl+MipQsxNeJ288SBR+Jyri7vlQTbKcggbfW1tJgxrCjdADzMf3C6zJUgmj/n58FpcddvsphGrVEOk/nrF+MFrMTQIpSxQxtsN79zvYmvhIqrZs8NpYdVcmPKD2W9jxRf8z19oMMwKUJtK0MBDWfrQmDgtcjpR7uJnmJ/vOqMxe8X4JcxCJjNq8Eu0Gfp3BfuaaFaNLtppZCrs3WkIsZ2oR5uNXiB+8orBuCG59YAXhgG64P/lDJUpDlJcv6EWRPJhFB1JlaonjSnsYVexTe10znePt+d/VWweVdzG1O/XxIe/bvTpcIvGm0cZvWEOPdNfMROE3v2xoRYWhWDYiaej0E4wctLDqDuxIlWqnjTm1njWznZh5Z7Pudq+reixN+SPu6jN/vk73zK+srrNY/3PvOeR2JG9AfyERi46EmqQaLO9eDp9NsxrPeyOhpX4fDcYjemxVwo/hUpefqO1x3TheH58ynh7tv+HiHp9EYK/XWlU1LBaVMT+W2O83bAM7ajG3gbDEL8ee2NG7NXELnqxYrVlrMf57W0NuAuEcv5gWshQug+3hxG0srdBu89UiuLFS7+bXYdKLO58beqMsfjo9K3TsaEYlYnMhZJpTI+9EvjNjZ/a6d1/Y+dmj4IOOxTMPM/9wzV6gvHBPOmOYg8yEoN4OxpFjx7Lz1NHluSVaw0tul9EvSS0pce5ClnUbQaLtxuL8YUet+N60RjjG6dT4BD24qQUGVUMidj94pYbXxip/PUxXKNqs6WuX6z9Nhx6e2Q0XLfnFKZt6GHedGcPcvaS1fT0JpgbWcO50Sf8zyGLUtBU65ZCLZNKZO20kLcznHlece8Wtwi68FtlJAYDu53nFcbipNWuiY8xzme2KzADJZDIPwSGr7AH9QSD9K1wE21rqwJH+8UlH9fG6slo+Hn6SH7WfxjxFS8BFIt71v7Z0juSuZE1LIndGWrEWBVboYZxfD+SFagajrcbmVIm5rW94kL7a9tu8XFgj7OouzXy8kthPYjduZwCJSEww1PYCy8yP9GC0+Km9eup+YFyimhVND/da8WF9se6feWd2gkb7obXnzGeqWFN5mCuSZ0HQCKZZumqTcyLPmv7Gj836iIuGOiiVQ04nqIUgoh6rd8E+nUUjaKFVFnzDOV1mjtgJfc4HKDy99H1e20XeuJOvjoqnE6Kfq1HafAVCsPTYy+l5YBX7myeB29DLJ5fsTe10/j6HLSA4+Gvw+vPQPaZSsEnIi9zbdNPcods6U04vsYrYssZp951Hr+KmJ0WMwX3unTIHRjdxtIaFiQvYptuK1uUy3p+JGZsS/etLcZ/jqsaOl+gS93L1ymYUJFwvG+/gi1FRqEwPIW9lHag5qKoWRSkosbPbrvBmzgVZpibUAcp4LDJdVcKzoquzv08ri3u+FpGJN7h3ab9nccfAtIOl12PHsvk/nu5NPll3s2MsrTUDXeh0W0sDazMzGac2mb/ew3b9SjPm4P/+Wb/7oWLnPNuz78O3IqArJTa6tZt4xe7JIGgC59+BFuKjEKjsYXd6eIrpST5xS6jf7WZRqjTxs/mmI4fHDXYVMxuPk7Nx5xwSGOMkjHPRk9vgncYa3vclsx+XJ84lT7tvMuQk2iFETVnNNyTPoJ+nV8126+j3JjqZG5kDYtid9GudqOUIZBRNXSLtTobGW/R9u/fDkbxkf5l/D9lf3MMJOpmFtNlf4KJh7sf67fa0ukaNjNenITYbJ1hV81cGPGX0q7abv7ZTK1Qq1IFoJGF3e3iK6Uk2esrrtfNws+HwU74Cx9T9n+yNJG8LOnvD8wnUSDefbqZG1OdPJiZzcLkBWzOjCWjFe9mRrFdjyKjFZszYx1zXkzB0hoSmWjJQr8+88+oAmtBoZgR+QtLYnfmRL0cUloxoFXgOUayr/7GVGfRza9PN7MoeS7xWJQ9B87B1kbzO29HG8/h2vC7LaTdtW3Xv2jFF+CGycXncOocag1cSu2eWjj/ebdneyT5DGoE3zTu4qnbxWcuBoWZFeNWVu01HzOaL1z06v6yoabWlDCH/jAPRo7KE+SVmdmQhCuauhin3mWL3o8bU53G43loEozgxuTg7/5rxMWMw96KAGNK7+oxPJmextnRJ4mgfQuxAm6K3UGTyheQZpXirOjqosdLRxOjdJ/b+v51qG2kiTCCARY2d3HJgdt5/5YHyV/4VnDomSi3rebMNEcVzRdCt77+bguedgTJeElsL15Y9ZMkUKrdI0VMQ0bjthQIu7Wnn/0X3VIoveYTpPVpbCSk9hq2jIry6sT5zNk01/dLMdMdrZkxWhs2w3fTn2NEU5Tv6duIuIhiRisuTX6JpbEf0qLs7aGgWSFDmS3jdK53M6OYMbAs97PTe2U7TzMtsKjyOFtcNvFwePAr+VvNRaKQcaoSDqkNreO1Z5m39Rq2C1Cs3w5kL9KqMbzb9kL4rT392DdufrnXfAKkeWWSfUxO/JxZIx7grHGPcGQAUQcjCi1Md1QK2tVuvhf9EXsGUvw8fWRRdoqVLXo/rmjqchX1u9JHGtaOz9gh4+JjhBV/aG0s0r7UPK1oTK3h4Uy+1+30Xtmyc7O95XDyMmOh/NEri/cPdRR1vK9VvwuYXuNYrz0/lo90V6x5GlfYK9G614/HWep8AtxwevVINMYi6bOvbg9cB+qU8QHQqga4oqmL9Zl/phdrZsogplfvNo4Grkmdx40pf++P1oP+ttN4KR3O5XpDspMJyVeLBFopOCG6Nu8xt9dYhPk3dLrBB602dbtWgyxgOmW8mKhIsAX9cj8LQsVpXI+9Eq09y/EIveZj9xU+EgMyRVHdSIyWAMV+uT+26LGMdxGsjsi2Ivsho7MZN3pszqu/Qnc5jmNmlVzR1OUY4VpvGF4WzBY9NpDImi23ClEKlsTuJK4HbA8Yw27AsGCuaOry3xMn7Ig13u5+rQVpcWH+7LQjl5lpFaSsv1y/XLpBVpTGjdgheCphNefjlDXQMrpomBaV5oomIzKbG1nDmuaLea3lTNY0X8zcyBrPadhlfFjJECmyHyLKEPXZA7flbig3pjqL0hYBBnRTLlJ3ywc30xm9RL1PN/NkZpr7QZZxnUTdxKvq9owRa1kSu5PxkW2+PH8N+TUNdgQp8jH7/LsRdAHT3KXr5B8NXmN+UhsrQSnpkkIgGlvY6w074Xfo5jhOvZtb2Bsf2UZEwfhspG2Ku5MmrcymO5pFQFb6VQvKIeWtsGJ1ZWY2C5IX5RUTbdej+Ebywpz4O+WD+10kTekIC5MXMCey0XUx1zpuWeuvCharH9iKv50tlX1Kfk2DHUHE0hRXt/FKXUOyXmN+UhsrQanVsYJvRNhrHYcPqrl4WShApkcORiTZ0RZ3FLoEI3LetdawOTOWBf3nO4rxFr1f0WMrM7OZMbCMS5JfpkePpY3dXNHUlbu53JjqZEDnO35+F0L7dDNfT37RtRK0VJxuLAqce9krGMBhO0I3YXqxK9hmz+AdxYaxhlStvUNLTZcUfCPCXuvMuZpUdETeQ7nFy4h9zxczsh7TGss9FrUoWWGk36QyJGjOeedOxTlOC6F23xxuid3O6y1nsih2Fwr3TT8KMSP/hckLPCP/Uig1w0ariGMWEDAoTNZslRsmGymOTqgoJdkiYSxgViu7RTajrjgi7LXO1E6uV1/MVYluzozNCd7fndoG6P2IRRW796bo6U2ggbRFzbwi/ZUFlanWc9phN15EDaZQxgpaAihV3Ogrk23+tTkzlkuSX+Yj/cvyzue1LmDipNkZyL2W7XqU5zhFxOJEnKwLk9Hji/3jxPbiFEfLmJx0R+m2SLlrSNXKbpF0yYrTuFkxDcT/3X0YP+OwoscXD8zn1pE/zfMrzcg6pbVjZOpka1g99JWZ2awc8Jd1U6pNsjkz1qMqdpD8SlrjfFbPXQM79CgeSh/O6dGn8iLrAaIsGLiIh/RsMtr4huFWWGWSImJsKmJmbbj10Lf22ffbb90qovEx9hkr8TH+xiqValSDymbUFUeEvQ4Y1xanp7dYLNbt+xn47HR48jr0zs1s0ftxQ9JdIME53dHOQ/eDV/qk03NmD9xGU0SRcquEsmC92ZwYWcOt+z8EOzfTFz+AX+36MEdENnBO9Al6GcXuDIxRexxvGoW9agqrSft0M1clL+DW7y/On4TdfrbxdiOLZWpntse+D0ZP8CdkiR2GpdNo4iftBSqKWDF1wIKjpxCP5Xuw8ViUBUdPyX0dnz1iBbP6b/OV2x7UQy9lPDe0Jpe+GIsWr2L6SeF8MDObWf23MWnvPVy18yROjT6T8/jb1W7iaoBLk1/KS880uaKpi2aVyntMKWMR2Wo9rdv3M/knta0q/ZGRRmjtteKFne3guJetRlIChaBIxF4HzJveAcDSVZvY0ptgXFucBUdPyT0O2c01fOKnQZhZoDNObWOLpSjJe7xim6QQpWBOZCPXAIlkvrdc2JtlvDJSOEmSd36zPTHAApc1Azs7yck6iqB5X/89gHHjXHz0lOKDvCJNp0Kzln0M8XaKvJ2ab1mRreMEn4iw1wnzpnfkCXkhTnaNE24eul9xdRpv8KbwLsqh86PTLk5uC7vm+Nb2xAAdDkLt9PjbjKXDpnvl39VYFNjeOH1Tqn9s2zzMBkkJFHwgwl7ndG/oYemqTfT0JooELxZVoCHp08M28SOublhFfk3zxa5+fqQgQ8ZrYbctHqM3kcz7XZoITRRnlWRUhHgsSiI5uEgaj0XZMuMKOl66pqiD4QEnfJ/Xpx7n+fo8KcU/LrwhqIh9Pr2kBAo+EI+9june0MNVK17KRerWUvqOtjhLTz2UpfMPzcth94OfrBm/uPn58Vi0KO3RqzhqZEsTHW35qXIRG1EHUDqDQjOmNYbCeE8Wn3wIH517UW02sbKmL550h6QECiUjEXsds3TVprxoFAarTZ9deETuscuWbww0bphZM65+fiZNVKm8HPsbU51FDcisC7tbehPcfNo0rlrxUu61O893LH3JDBrFzadNy7dWaj0rQ1IChTIoS9iVUvOBRcAHgcO01hXePUOw4rRgWvh4UP/dS1yD4ubnpwuS7b0Wdse1xYsWk+9sPpsrkrc7zjeRTLN01abSPPNqUus3H6FmKTdi/xNwMvDDEOYiBMRJsNssrQTASJe0Rrhe+N9WrzI43QhyKZ4ULiYfx6Lr4YKBux3nGyRrSBDqnbKEXWv9ZwA1VPuZDXPMhVIz5fHTB+3PfX94i3SBUb17b4ruDT054SuMcNtaY2gNOxPJ3DhPvbI17yYRpPK0kphWTYdHpsq04y7kMys+RmLA/uY1rs1lowlBaDDEY68TzIVSM+ru6U2w3EbUwciCKbQevNIluzf0cGlALz5MCjN64rEoi08+xLd9Yh63aOXLRVkz1khfEIYDnptZK6WeAA6w+dW3tNYPZo95GviGm8eulLoQuBBg4sSJM954441S59yQFEbjhdHprCWrA/nkCnh9iX3qnt25zJTJahCPRTllRgdPvbKVLb0JRsdjJNMZ9mSj77Z4jEVzD/Yt8l7vpSDUK343s/YUdp8nexoPYbcyc+ZMvW6drLOaFEbjUByxTl74SOC9Te3sC6dz+fXfw6Zwjt0beljwqxdIpvNfbSyiWDr/0CKBFhEXhhN+hV3y2GsAu7RFM5PDpBSPuKc3wVUrXqJ7Q4/nuaqxStIWNxZ5L1u+kVlLVudEulDUwbCXLrUcB/l5/Obm3oWvVxCGI2UJu1LqJKXUZuBjwCNKqVXhTGt44Sdt0a4RWCyibJtoWSm8QTidq/zvbcGIRRR7BlJ5onzp8o2edlBPb4IFv3yB6dc9xqXLN3reEAVhOFJuVswDwAMhzWXY4pS2aI3SnRqBWR9zEmermLe1xtjRl3Q4srKYC6QdbXH6BlIlzyOZ0a7PldRGYbgjWTE1gF2euV0mh1Nmi/mY0wJrW2uMWUtWlyR4zVHFgI01EoSzD5/I9fMOyXts8sJHyhrTDUltFIY7IuwVxs/inp+2vH7O4dQIbPfe0qLjWESFYtHcv76HmQe2572eoNWwfpHURkEIKSsmKMMlK8ZPtkslzmG1PPb0p4ryuv1Qrl1iN561f43dvMM4h2TFCI2MZMXUAH6yXZzo3tDDrCWrmbzwkbxMkEIWrXzZtRHYzhJEPRZRLDh6Cr0hevGFNtC86R0sPvkQxhS0P7CdT1TRFjc6NLbFY0ULxvFYlFtOm8azC48QURcExIqpKH6bdBViV2V62fKNXLp8Y15U2r2hxzEaN8/htVhaaN3AYOWq34VWs+y/sFOjFTvf21wzKLSrJu0XZ+1rO3JjnvbRCXkeveSuC4I7IuwVxE+2ix1O7XhhMFfbPM7t3N0beti9N+V4jFthUk9vgpjb/naWMazWkpP95OZ7WxeFzeebN4i01ty99k0efuHtXPWpV3sEQRjuiBVTQVw3oXbBK6I37Ry348w2AW67JyWSadtt60ycnmtu3GFuXFHYk2bxyYfQ0RbP29zCrxDb3dQAehNJKT4SBJ9IxF5BSs128ZMxYo5nd9zI5ijzpnf42mAj6Nq5Al5d/Flg0BK5bPnGvNdWTkTtdrOq277qgjDEiLBXmFJEzk//dFNI7fqqDKQydG/oqUhKoWkj2a0DmBZROf1cvOYsxUeC4I1YMTWI1c4Aivq4KODTB+3PvOkdtj64ufhpZwWVg9VG8pvxE7Sfi9ecpfhIELwRYa9R5k3v4NmFR/A/S47jrMMn5om7xij6OetHz9GXtN/IeUtvIneDCINCr9xvxk/QlM950zs4ZUaHrfcvxUeC4A+xYuqAR158uyglMZFM8+yr2x2fMzo+2EbALQ3RD4XFReA/4ydoymf3hh7uX99T5P2PaY1xzQn+e7ILwnBGIvYK4Le4yO9YpVR/WjsnliPqYC/CfjN+nKwTp8edsmJam5tE1AXBJyLsIRN2j/BSWtAqhW1PczfisWiuP3ohTsVFftIag6Z8llrUJQjCICLsIVNOGwE7ggqaIngKI0BLU4TjD31vSXn3bgTNaw8a4QuCUIx47CETdsTp5GXHYxFAFTX/OuvwiTz1ytbAaY69iSTL//AWzU2D93q3vUaDpDsGSfn028JYEARnJGIPmTAjzu4NPfQNFLcEMMr4pxZFwjefNo3r5x1ia3/42foumdG5DaQB+lP2GTcA1z5U3HwsjN2Lyq1cFQRBIvbQCSvidGprWxhFu228YS0K+vRB+3P/+p5AbXKdKj3dFnTD8MKlF4wglIcIe8iUu2mGiVN2yMiW/OwQp6pOO3GceWB77tiIzxRIO6H2aj4mCEJ1EWGvAGFEnH68+m93v8Q9a9/M6/xo197Xbl5+N7qwE2qv5mOCIFQX8dhrFC+vvntDT56omxS293VKsyz0su02sADoG0gVjeE0t7Z4TCwUQagBRNhrFK/876WrNnnuR+q1mGm2LXh9yXFsvOYolp56aFEu+46+4na5C46eUtSjJhZRLJp7sI9XJghCpRFhr1G8skP8pjMGWcycN72DkS3F7pztDcKuM5kgCDWBeOxVxKudrZNX372hx3ZLOzuCLmb68faXrtpUVNmaTGvplS4INYJE7FWinNYDfmwYKC3N0k8evpT9C0JtI8JeJbxaD7g1EnMT0HILe/z0dpGyf0GobcSKqRJuUa9Xub5TmwG79rpe2NlBi08+xNUikrJ/QahtRNirhFs/c7doft70jopVt5o3kMUnH+J6gwirCEsQhMogwl4l3MTZaRNqM8qvZHWr3w2jpexfEGoXEfYq4SbOS1dt8tydaKiqWwVBqD9k8bSKmLbKuLY4W3oTLF21ie4NPYE3pygVWQQVhMZEhL2KOKU8AkPSunaobiCCIAwtYsVUETeP+9mFR5Qk5F5FT1ZkEVQQGpOyhF0ptRQ4ARgAXgU+r7XuDWNiw4GwPe4guxqZyCKoIDQe5UbsjwNXaa1TSqkbgKuAK8uf1vDALeXRL9YI3a7Hut8sF0EQGoeyPHat9WNaa3PvtrXA+PKnNHwo1+Mu9OidNs6QLBdBGF6E6bGfBywPcbyGp1yP22mXpUIky0UQhheewq6UegI4wOZX39JaP5g95ltACrjHZZwLgQsBJk6cWNJkG40gC512+InEJctFEIYfnsKutT7S7fdKqc8BxwNztHbeRFNrvQxYBjBz5kw/zQkbmlIWOgtx8uijSpHRWrJcBGGYUm5WzDEYi6Wf1Fr3hTOl4UE55fwmTm0JKpHzLghC/VCux/6fQAvwuFIKYK3W+otlz2oYEEaqo+ShC4JgR1nCrrX+p7AmMtwII9URJA9dEIRipKVAlRjqcn63jTsEQWgspKVAlRhKGyWMhVpBEOoHEfYqMlQ2ShgLtYIg1A9ixQwDpO+6IAwvRNiHAdJ3XRCGFyLswwDpuy4Iwwvx2IcBku8uCMMLEfZhguS7C8LwQawYQRCEBkOEXRAEocEQYRcEQWgwRNgFQRAaDBF2QRCEBqMhs2LK3ZlIEAShnmk4YZeGV4IgDHcazopxa3glCIIwHGg4YZeGV4IgDHcaTtil4ZUgCMOdhhN2aXglCMJwp+EWT6XhlSAIw52GE3aQhleCIAxvGs6KEQRBGO6IsAuCIDQYIuyCIAgNhgi7IAhCgyHCLgiC0GCIsAuCIDQYSms99CdVaivwRpnDjAW2hTCdoaTe5lxv84X6m3O9zRfqb871Nl9wnvOBWuv9vZ5cFWEPA6XUOq31zGrPIwj1Nud6my/U35zrbb5Qf3Out/lC+XMWK0YQBKHBEGEXBEFoMOpZ2JdVewIlUG9zrrf5Qv3Nud7mC/U353qbL5Q557r12AVBEAR76jliFwRBEGyoa2FXSn1XKfWiUmqjUuoxpdS4as/JDaXUUqXUK9k5P6CUaqv2nLxQSs1XSr2slMoopWo2s0ApdYxSapNS6m9KqYXVno8XSqmfKKX+rpT6U7Xn4gel1ASl1FNKqT9nr4dLqj0nL5RSI5RSf1BKvZCd87XVnpMflFJRpdQGpdTDpY5R18IOLNVaT9VaTwMeBq6u9oQ8eBz4sNZ6KvAX4Koqz8cPfwJOBn5b7Yk4oZSKAj8AjgU+BJyhlPpQdWflyc+AY6o9iQCkgMu11h8EDge+UgfvcT9whNb6UGAacIxS6vAqz8kPlwB/LmeAuhZ2rfU/LD+OBGp6wUBr/ZjWOpX9cS0wvprz8YPW+s9a61rfCfww4G9a69e01gPAL4ATqzwnV7TWvwW2V3seftFav621/mP237swhKemNz3QBruzP8ay/9W0RiilxgPHAXeWM05dCzuAUup7Sqm3gLOo/YjdynnAo9WeRIPQAbxl+XkzNS469YxSahIwHfh9dWfiTdbW2Aj8HXhca5vAH30AAAHUSURBVF3rc74FuALIlDNIzQu7UuoJpdSfbP47EUBr/S2t9QTgHuCr1Z2t93yzx3wL46vtPdWb6SB+5lzjKJvHajoyq1eUUqOA+4FLC74x1yRa63TWqh0PHKaU+nC15+SEUup44O9a6/XljlXzW+NprY/0eei9wCPANRWcjide81VKfQ44HpijayTXNMB7XKtsBiZYfh4PbKnSXBoWpVQMQ9Tv0VqvqPZ8gqC17lVKPY2xrlGrC9azgLlKqc8CI4B9lVJ3a63PDjpQzUfsbiilPmD5cS7wSrXm4gel1DHAlcBcrXVftefTQDwPfEApNVkp1QycDqys8pwaCqWUAn4M/FlrfVO15+MHpdT+ZuaZUioOHEkNa4TW+iqt9Xit9SSMa3h1KaIOdS7swJKsZfAicBTGanIt85/APsDj2RTNO6o9IS+UUicppTYDHwMeUUqtqvacCskuSH8VWIWxqNeltX65urNyRyl1H/AcMEUptVkpdX615+TBLOAc4IjstbsxG1nWMu8Fnsrqw/MYHnvJKYT1hFSeCoIgNBj1HrELgiAIBYiwC4IgNBgi7IIgCA2GCLsgCEKDIcIuCILQYIiwC4IgNBgi7IIgCA2GCLsgCEKD8f8BB9gIdCPFWhsAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(data[:,0], data[:,1]) # Print the original data in blue\n",
    "\n",
    "# Apply PCA. In theory, the Eigenvector matrix must be the \n",
    "# inverse of the original rotationMatrix. \n",
    "pca = PCA(n_components=2)  # Instantiate a PCA. Choose to get 2 output variables\n",
    "\n",
    "# Create the transformation model for this data. Internally it gets the rotation \n",
    "# matrix and the explained variance\n",
    "pcaTr = pca.fit(data)\n",
    "\n",
    "# Create an array with the transformed data\n",
    "dataPCA = pcaTr.transform(data)\n",
    "\n",
    "print('Eigenvectors or principal component: First row must be in the direction of [1, n]')\n",
    "print(pcaTr.components_)\n",
    "\n",
    "print()\n",
    "print('Eigenvalues or explained variance')\n",
    "print(pcaTr.explained_variance_)\n",
    "\n",
    "# Print the rotated data\n",
    "plt.scatter(dataPCA[:,0], dataPCA[:,1])\n",
    "\n",
    "# Plot the first component axe. Use the explained variance to scale the vector\n",
    "plt.plot([0, rotationMatrix[0][0] * std1 * 3], [0, rotationMatrix[0][1] * std1 * 3], 'k-', color='red')\n",
    "# Plot the second component axe. Use the explained variance to scale the vector\n",
    "plt.plot([0, rotationMatrix[1][0] * std2 * 3], [0, rotationMatrix[1][1] * std2 * 3], 'k-', color='green')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The explanation of this chart is as follows:\n",
    "* The rotation matrix used to create our correlated variables took the original uncorrelated variables `x` and `y` and transformed them into the blue points.\n",
    "* The PCA transformation finds out the rotation matrix used to create our correlated variables (blue points). Using the PCA model to transform our data, puts back the variables as our original uncorrelated variables.\n",
    "* The explained Variance of the PCA is \n",
    "\n",
    "$$[1.0094, 0.1125] $$\n",
    "\n",
    "which is approximately\n",
    "\n",
    "$$[1, 0.333 * 0.333] = [std1^2, std2^2],$$\n",
    "\n",
    "the parameters of our original random variables x and y\n",
    "\n",
    "You can use the previous code to try with other standard deviations and correlations and convince your self of this fact.   "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PCA as a strategy for dimensionality reduction\n",
    "\n",
    "The principal components contained in the rotation matrix, are decreasingly sorted depending on its explained Variance. It usually means that the first components retain most of the power of the data to explain the patterns that **generalize** the data. Nevertheless, for some applications, we are interested in the patterns that explain much less Variance, for example, in novelty detection. \n",
    "\n",
    "In the next figure, we can see the original data and its corresponding projection over the first and second principal components. In other words, data comprised of a single variable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnX+QHOWd3p/vzPZKs0uOkWJVYQYWccQl5WSBhDa2YuXqAnaQfRh5LdnIDiSu3FVRrlwqElFtLMCHVhwulFL5QMldKuFil+sKHZZAeE9YvhK2RcUOnChLrBZZRnIMGMEAOZ2lhUM7QqPZN3/M9mxPz/v2j5me6Z6e51O1pd2enrffntY8/fb3pyilQAghJD1k4p4AIYSQaKGwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyuiL46Af+tCH1OLFi+M4NCGEdC1Hjx79e6XUIr/9YhH2xYsX48iRI3EcmhBCuhYReT3IfjTFEEJIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyqCwE0JIyoglQYmQNHHg1QPY9eIuvHP+HVwxeAU23bgJt/72rXFPi0TI+EQROw+ewltTJVyZz2F07RKMrCzEPS0jFHZCWuDAqwcw9vwYLlQuAADePv82xp4fAwCKe0oYnyjinqeOo1SuAACKUyXc89RxAEisuNMUQ0gL7HpxV03UbS5ULmDXi7timhGJmp0HT9VE3aZUrmDnwVMxzcgfCjshLfDO+XdCbSfdx1tTpVDbkwCFnZAWuGLwilDbSfdxZT4XansSoLAT0gKbbtyE+dn5ddvmZ+dj042bYpoRiZrRtUuQs7J123JWFqNrl8Q0I3/oPCWkBWwHKaNi0ovtIO2mqBhRSrU2gMjVAP4SwBUAZgA8qpTy9BwNDw8r1mMnhJBwiMhRpdSw335RrNgvAdiilHpRRP4RgKMi8kOl1C8iGJsQQkhIWraxK6XeVkq9OPv7PwB4GUByn1EIISTlRGpjF5HFAFYCeEHz2l0A7gKAoaGhKA9LSKww8zRaui3LM4lEJuwichmAfQA2K6Xec7+ulHoUwKNA1cYe1XEJiZMDrx7AvT+9FzOYAVDNPL33p/cC6O7MU1tci1MlZEVQUQqFDohsN2Z5JpFIwh1FxEJV1HcrpZ6KYkxCuoHtz2+vibrNDGaw/fntMc2odWxxLc4m4FRmAyxskR2fKLbt2FFneY5PFLFmxyFcu/UA1uw41Na5J4mWV+wiIgC+BeBlpdSftj4lQrqHUkWffWja3g3oxNXGFtlWV88mc0uUWZ5Rrv67zTwUhSlmDYB/A+C4iByb3XavUuoHEYxNCOkwfiJanCphzY5DTYucTnDv3nMMR14/iyvzudqTgpNmsjy9Vv+tzjfp5qEoomL+j1JKlFLXK6VWzP5Q1AnpUvxEVFAVN4XmzDM6wVUAdh8+jZuWLoosyzOq1T+LgBFCuh5dCr2NoCrCTsKKnElYFYADL72N+dacLOVzFh5av7yplbHpBpURCXUjYhEwQkjXM7KygIfWL0dhVhizIgCAQj7XIOo2TpHzc1h6PRGcmy7j3HS59vcHl2aM+/phukFVlAr1lNGNRcBYK4aQLqJTTryRlQXtuGt2HPK0gXvZowHUQiiD0oqz1n7Plr2TtcieZsYdXbuk7pyA5BcB44qdkC7BGYbYrH27VfwqHZrs0WP7T9SFUIahFZPHyMoCZgz1sIKO63yCEVSfXJo1D3UKrtgJ6RKiivJoBb9KhyaxnCqVtduzIphRClfmczj/wSXtfq2aPKKItDE9wSQVrtgJ6RKS4MTzMwWFFeGKUnh44wo8t/VmjK1b1pa6591YT71VKOyEdAlxO/GCmIJMIrpgwDKOa4/RLpNHN5pSWoWmGEK6hLBOvKgdrUFMQSZTDYCGuevGiNLkkYRs0bjmQGEnpEsI08mnHdmSJpNPcaqEFdufgQgwNV02zuvI62fx2OHTocZuliRki8Y5Bwo7IV1E0BVtOxytJickUO8cNQnYsyfPeI4dJe10NAddhcfp7KaNnZAU0g5Hq1dGqhtdNqrXsaN2ZLbL0Rwm5DROZzeFnZAUEtTRGqasrTsj1Y/iVKluPNOc8jkr8hWs8VgeTtwghKkbE6ezm8JOSBfiJ8hBQvyaSXgaWVnAc1tvDizuo09MYuUDz+DarQcwffESrIw0zGls3bJAY4VhdO0SWFlp2P7+hUstJXSFWYXHGWZJGzshXUYQp5xXdIpdcjcz2xnJSalcwZa9k3Vj6NBF6Ogoz6ha7Rf735yVwYXyTGRRIiab99j+Ew0JT+UZ1RY/g24VHsbZHTUUdkK6jKBOObej1X1DcIu6jV0kyx5Dh719855j2te9uFCewcMbV0QicF43uXcN2a6t+hnChJzGlbFKUwwhXUazTjmvzkhudHZjt/kHQGCTjBMFYGz/idDv0+F1k2uHjbtbkp24YiekyzCZA/IDlmdno7ArVdv5ObKyYOx6pKCv0e7HVKlcG7sVvG5yD29c4bm6bjZ5qBvqxnDFTkiXoXPKWVnB+xcueTpCm1mp2mOYuh45/3WTs7zlJYoORF6rcq/VdRIqZbYTCjshXYZOsAb7+1Ce0dcctzFFady5esgYn26P0YxdeuHgPNy5esj4ujscshn8Ik/sKJ7XdtyK57beXOfQ7LZ2d2GgsBPShdiC9fDGFQDMZXGdgmzfEPK5uVju+VYGw9csxEPrlxuPZZsqwvLWVAkPjiwPVACsFeb1zcnYgoFgrfS8yiMEjetPMhR2QroUpznBhE6Qne3mzk2Xa1EkJgG27c+NUeHe2MfedltjOV6bVlbJ9vk7b2oXysFa6ZluVK026k4KFHZCuhS/KBddGJ5Xh6P3L1xqGMPKSs2pGMZB6jaH+D0RNEMr5hSdCSeKRt1JgcJOSJfiJYimMDyvDkduGz0ADPb31cYwregzAty5eqjO5r9hVQE7D54KFBrZbPihlznFz5Si81MEadTdLTDckZCEYwrLM4U9FvI5PLf1Zu04umxTL6ZKZazZcQija5fA9Lbfmm/hwZG5FbkuNHLznmMYsDKwMtJwA3nr3RIWbz2AQsjMTK9qk05TCqBPtHKHLfo16u4muGInpMOEKbzlFZZnigi5aemihvHtccKIuo19TJOD1hZ/+zxMJqLp8gwgjWGQ9pTC2rRH1y5pqD3jxmlKiaK+TrfAFTshHcSvzot7dT598ZLRjmyvyp3737R0EfYdLTaMP9/KBM461VEqV5D1WO07z8PLdFGuKMx4+DfD1CsfWVnA9qdP1GrQmHhrNqyy2fo6SU9G0kFhJ6SD+Dn83OJjwn5NZ07QjR9E1PuzgosV84q+ohRyVtY4ljOV32vufk8NYWzaUz6iDlRNKc3W1+lWaIohpIN4pcCHqeUiQKjmDkHwEnVgziHrVR/mramSb0OOrHibT8LYtP32tU0pcTa9iAMKOyEdxCsFPozIKOhT8r2aWQTtfmTCNkt41WN3pvIP9uuPlxGlrZUO+DfndtvITWGLQH1kUJxNL+KAwk5IB/Fy0IUVGV1Yn2n8sXXL8ND65b6rZS+cJoogqfz5gX7tOOUZAGoufNKeU87K4INLFWzecwzX3fMDfH38eO094xNFjD4xWedEHn2iWjfeHbb48MYV+LWrhECaHKNBoI2dkAA0WwnQjZ+DTleNcF5fxhiRYgrr85prMzXU3TeEIMfxdKLOhjwWHPb4kiNrtKIUHjt8GgDw4MhyjO0/0RAmWZ5RGNt/Ase23dJwLXTX66H1y1PhGA0ChZ0QH7wiKprB5KCztzkjPeb1ZfDZGz5cF+miw935yKtBRpBIEjdf/vjVgc/Dxs+Jem667DuPx194Aw+OLPcMtXRjul4PrV+uje9PI5EIu4h8G8BnAfydUuqjUYxJSFLwjGT5cPTHc9Y7mSqVse9oERtWFfDsyTO+0SbOG45pdbrttmXatnamuuqD/dm6BCQ37tXxTUsX1ebaTK129zmFJWgETJqJasX+HQB/BuAvIxqPkMTgFVFxWQTC7hRGUx/Sx194A9+8/Qbf1bZd9+X8xUsoz0a5FKdKGH1ysraPLXxu0dVJaM7K4huf9xb10Scma2aS4lSpZkIxjRkG2wS0YMDSnreuzEGvRcDoiETYlVI/EZHFUYxFSNLwamD8Xotjh+lDGtQ2rjNPlCsK9z71EhSkdiw/0XWm+Hs1jNbVmHGzYMDCQH9f7f3nP7hkNK84sU1A225bhtEnJ2s3K6BaoGzbbcsa3hOm4XRa6ZiNXUTuAnAXAAwNmYvvE9Jp/ByjXg2M//il1o4dJna9VaYDlrQF6uvNePkYgogzUE0kmrj/ltrf7jHdZAT41x8fqpmAwmSFhm04nUY6JuxKqUcBPAoAw8PDrT6hERIJraaatyrsXjbzOHGaLbxK/QbFuVp2ttqzyxTY/3oVAguaFZqm0gDNwqgY0tOYRGvznmPYefBUTRDalWruVX8lTpxCbLr5BF2t26vl8YkixvafqHufXaYgSNcjHaanrbSUBmgWCjvpabwcan5lX71QCoFK0SZR1K2s4PwHl3Dt1gO4Mp9DRoAAZnQtAmDDKn2Mvk2zEStBnrZ6lUgyT0XkcQB/C2CJiLwpIn8YxbiEtBs/h1qrHXS8StGOTxRbygRtG6q6GreTn5oV9dmh8OzJM76+hGYiVtLekLoVooqK+XIU4xDSaXSONjethsnpVqRe9dG9KihGwYIBqyrcGsEWIFCUSxiCfH5BIlbcZheTiaiXwhpN0BRDehqno80kFFGEyTnFZnyiiC17J7WinhXBjUOX47lXzrZ8TBPnpsvIiD7cMaikZ6Q61yA3AfvzM32+dnOQlQ88U4tVz+csjK1bVrs+OrOLKfmpl8IaTbAIGOl57IqFj2xc0bZCUbbY+HUyqijVVlG3aXVRfnnOws4v3lArvrVgwIKumZGVqTbDNpXyXTBgYcOqAvb87I26BKSpUhmjT0x6dmVSmKvkaNNrYY0muGInZJZ2hck5xaaTcevt5Nx0ue5zOv/BJe3N4rL5fXWfn+6zXbPjUF3ikU15RtVMWCbzigIChUr2GhR2QhxEHSa3YMDCttvmTAppsv/aphWvWHxnhyPTZ+v1mdivednU7ZBJivocNMUQEhK74cPirQe0DkgnU9NlbN5zrFYzvRX7r64uStIJcr5e+9iv+XVlYjRMPRR2QuDfwd653z1PHQ+cMWrrvh32eNPSRdqOPznL+6u45rqFmLj/lgabcpJxJiZ5fbaja5doOyrZ9nkAta5Mfm35SBUKO+l5nGLtbFyhE/dWbOSlcgXfn3wbG1YV6gRaAbjk48389W/mTBJJZcGAVdfJ6KH11Tovfp/tyMoCdn7hhronkvysc9bdaNqvLR+pQhs76XnC1O9udVU4VSrXlbW1KVcUBvuzOH9Rf9Ow2+DlByxkAAQv59UZcla2zpdgs2bHoUCfbRjfBot8+UNhJ4klbDu6ZtvXhanf7dcVqBXOX2yske5EAaE7H3UCr2iUqGqju6+t3XikV4t8+UFhJ4kkbB2QZvb3am4B6B/tg2SqtkLyKsf4c9PSRdh58BTu3nOsQWSjqI2uu7b7jhabLhzWC9DGThJJ2DogYfZ329RNaf32o73T+bfz4ClsWFWo2ZK7yZnZLh47fNpoQ9dFs4Q1m7AmTHgo7CSRhH2ED7N9+9MnPFfc+ZxVWw3qHKv7jhYxunYJXttxK5Vdg1N0ndEsTqdqmJU2W92Fh6YYkkjCPsJfnrO09cHdYYTjE0VfO/XgvL66LNReb4zcDE7RbTXpi63uwsMVO0kkYR/hTdVvp8szWOyInw7y+O4UJdOqsDhVqpobutEo3gGiFN0ozDm9BlfsJJGErdvitwq3bb9BnJ5OUTI9CQDA6JOTmPcRWmPcuP0TrdbeYau78FDYSWIJ+ghvyhJ1UypXIALPMgDulaBXH4xyRWFeoCMnkztXD+HxF96ItItTVqTBPxFFh6Neb3UXFgo7iYUoVnI2258O3lRZqWrrN2c1QTt2vJDP1UL3Nu85lth+pFGx72gRX/741dh9+HQkFiUrK9j5hRvon0gAFHbScaLuVRk2aWewvw+D8/oabirjE0WMPjFZax6RZlEHqiL7+Atv4BPXLcTzr5xtSdzdVSwBc9XHdiV4kTko7KTjxL2Se7dUxrFttwCYe3K4e88xAL3nC60ohRdPv4tPXLcQh189F/pmJgDuWD2EB0eWN7xmeuJJZJ/XlEFhJx0nSFxyGFNN3uDg9Gud5n5y6FVK5UrTK3YFYPfh0xi+ZmHD9fHqEgVEa44j9TDckXQcUyicW3CDVFsEgLF1y2C5+rJZGcEdq4c8w+TS0s0oClp5UlGANozUVIWxkM+FvsYkHBR20nH84pLDppCPrCzU9d8s5HPY+cUb8ODIcm3WI1CtOkhbbzh0PU1tdE9hXteZZQLaC00xpOOY4pIBb8H1SiE3hcO5t9P80hwLBixPJ/WVs6tw9zV9aP1yrbnF9mm4YZmAaKCwk1hoRnBbzWYcnyhiy97J1Ee7hMEuZOZV392utb7z4CntTVdQrfCoi3R6aP1yPLf15ob3sExAe6EphiQCP3t3qynk9o2Dol6PgreoO4t26UwrdlTMsyfPhDKtsExAe+GKnSQCr0dwr0YOQaGjtDmcq22v1P5rtx7Qvt90XVkmoL1Q2ElbCBvKZno0L+Rz2kf5MMcYnyjSUdoki2cFO5+zMLZumdGX0YxphWUC2gdNMSRymgllC/toHvQY9n6kNaZKZYw+MWm8hjStJAuu2EnkNJNZGvbR3HSM7U+fqBtj+uIlmmAiojyjjNcw7PVjclJ7obCTyGm2402YR3PTWOemy7WwPJpfoqeZkFM3UdcKIo3QFEMixy+ztJ3HIO0lis+dyUnth8JOIqcT9lbdMUjz5HNW7ffB/qw2y9TKSCTXkD1M24+oCOJ6ReTTAHYByAL4X0qpHV77Dw8PqyNHjjR/wJf2Aj9+AHj3TeDyq4BP3g9cf3vz47UL9zw/cgvwf58xz/ulvcDffA0onZ3bZg0CffOA0rm595w+DBz9DqAC2I4lG2y/iFHOX2ZFotWafroxG/7W/Xc2VQOLgOuvvVrfjUMpvPTaG+05aNRI47VR7l90OzR5XVXdAbzngdxC4DP/pfr7+B8BMxfnXsv01/9t0z8IXDwPzwufWwgs+3zj9xHQa8v3/9Pcd04yQF8OKE/760/EWiUiR5VSw777tSrsIpIF8EsA/wrAmwB+BuDLSqlfmN7TkrC/tBd4+j8CZcfd3coBt/3XZIm7bp5unPN+aS8w/u+BGZ/a4jEJNdGzfLFZ2I//ukuEPfFk4J1GFRHZ/monFud30MoBV30MeO1/m99n0p82aFVQYY/CFPMxAL9SSr2qlLoI4LsAPhfBuHp+/ECjWJZL1e1JQjdPN855//gBf1EHKOqkB+mAqANA5WLjd7Bc8hZ1ex+d/sSoVVEIewGAc2ny5uy2OkTkLhE5IiJHzpw50/zR3n0z3Pa4CDofe7+kzZ8QEhzd9zdGrYpC2HUmtgb7jlLqUaXUsFJqeNGiRc0f7fKrwm2Pi6DzsfdL2vwJIcHRfX9j1KoohP1NAFc7/r4KwFsRjKvnk/dX7VROrNyc4yMp6ObpxjnvT94PZCzv/YGqjZ2QnqJDwXvZ/sbvoJUDrv097/eZ9CdGrYriE/sZgI+IyLUi0g/gSwD2RzCunutvrzofLr8agFT/TZrjFNDPc/gPzfO+/nZg5L9XvfVOrMHZbbPv+fz/qI4TVOA7dCNQqPqdGn5M24P+eI2dgB/Pz6QTc3BfA7/jBt2nieO3hdxCYP3/BNb/RTUKxon7b5v+wdlfPOJ1cgsbv4+f+/Pqd9D9Hf3K/vrvnGSq30s//YlRq6IKd/x9AI+gGu74baXUN7z2bznckSQSXZr43XuOtfTlt4uAfX38OB47fDqyuUbFZUu3moJi8P5Jz6jf0DiD9+yiXACMddLdWFnBYH+ftj8sUG0y/c3bbwg0Xs7K1sr5ks4RNComkpICSqkfAPhBFGOR7sWZUm6LfKvLhuJUCXf8xd/iuVfO+u+ccuzPsuDoOBW0G1RGgEpFGUUdAGaU8uxuZJMVwYZVrMyYZFgrhkRO1O3nKOr1FKdK2LznGDICzAS8cwbZzy4XYCrBa1NRCvuOFjF8zUKKe0JhSQESOV5NLRYMWBiw+N8uCoKKehCcJR+ClGtgbZdkwxU7iRxTzQ8BMHH/LQCAr48fx+7Dp9vvfCO+uDtUuUvwmq4Ra7skFwo7iZwg3XQeHFmO4WsW1sQjZ2UwXa7PMGxjiRcCbweo01+yZschNp7uMvhMTCJndO0SWNnGUJHzH1yq68AzsrKA57bejNd23Ipf/Mln8MjGFXVVBpMi6gJgzXULYelKHiaIASuDfM6CoBo1474GVlZqrzubVPvB7kjdB1fsJDC6cEagsWsOAK0qT5XK2LznGP7zk5MY6O/Du6VyXfecI6+f9YzaiAsF4PlXzibmRmOiVJ7BdHmmLmomii5FbDzdfUQSxx4WxrF3nlZbkekiXaysAKraMs0mZ2Ux38rUuhgFoZtNLp2MYw8D48zTSUfj2EmyiaIVmS7SpVxplONSuRI6zLFbRT3J+PWYBdh3NM3Qxt4DRNGKjBEQ7aUd3aC8rpl9sy/ORr3YN3unD4R0LxT2HiCKVmSMgGgfWRHcOHR55ON6XTP2HU03NMX0AEHCD924H9NvWroI+44WI8smJXNUlIo8u9YvaoV9R9MNV+w9QNhwNd1j+r6jRWxYVcCCgQClhR3kcxZyLWSaJjvAMLn4OU5NN3U+maUDCnsPMLKygIfWL0chn6vFMG9YVcDOg6dw7dYDWLPjUJ1t1fSY/uzJMxjoD/6Ql7OyGFu3DC//yWcw2N+cDZmO1fAU8jlfJyhj09MNTTE9grvyoleUTKuP6QLUxbSv2XEI5y/ShNMJgoozY9PTDYW9B/FynI2sLBht8hkRVHzyHhYMWLV6MFFXeUwrOSuDUrn5hs12HkBWpM4B6ifSzps9SRc0xaSQ8Yki1uw4pDWzjE8UjSVZ35oqYXyiiOmLl7Sv+4k6ALx/Ya5sgFeVR1Klmkh0Pe5cPRTan2Cb1e5YPYScla1dH4YuEq7YU4aXmeXI62ex26MLUX7A0q6ww2SGlmdUbeXfCxEWrTh3nVUVR1YW8OzJM4E6Idm8tuNWAFVTl+4JbMveSQDBk9BIeqCwpwyTmWX70ycwNV02CrSVFSgF7Qo7rAPTFvT8gBWqtEA3otC8uNsNM8b2n8DYumWhboQFR/SK6X0VpUJnGJN0QFNMyjB9yc95iDoAQCGyAlxX5nMYnyji/Qt6k04vIIK6SpVe2MXRTNfHfeNwO0i9QhSZdNSbUNhTRrNxyOUZhayumpUB056Caijd9qdP1BUH60WObbulKdu5k5yVxR2rh+pCVd0x6n4dj3rBJEbqoSkmZYyuXdJgJ89ZWczry/iuyCtKIWdlAzk8P3HdQm225CeuWwgAqTfBBMXZUCSM/Rxo7Gxkwn59y95JrYObSUe9B1fsKcGOhLl7zzHM68tgwUB9Q4Wxdct8C03ZiUt+K/dCPodf/0YvUodfPYfNPl3uew27oUihCYG9e8+xhsgm0zG+efsNTDoiALhiTwXuSJipUhk5K4uHN65oWO3ZK0d3pEvOytbqwXiFNdpCcbdBvIOERHYr/VlBZab5cxxduwSjT05qyx27EaC2wg9aZplJR8SGwp4C/BKObNzZp24BMMWdZwSwzeXz+qoPeaYkpjRz+z+7Gg+OLAcw9/m9a9pZVfcZWVmo+6zzAxY+KFca+rs60YWXBqmvDjDpiFShKSYFRFWpz7S/0wc6VSrjnqeO46aliyKvIR40iiQsLdQgq+PZk2dqv9vmFZPVSgG456nj+Pr48bqCatXoJG9Tl2k9TycoCQqFPQWErdRnarJweUBhtQuCOQuLeRE02mZwXl/o6pFB2PixIQxEoO52Zq4zq9eLUrmCx194Q/s0ZfpMCvmc0RZPJygJCoU9BYSt1Gcy3YgET7axzTDPbb25lgFp4rdywSx+xalSW2Lfnz15xliLRQCsmY3k8cPOzHXeEP3M7SZ7vB2B5MS+Zqy8SFqFNvYUENZpZnqk98pM1TH6xCSOvH62zkThRhAu9LEdse/2Z6LzCVyes/DF4SEce+NdzwqUOStrzMz1ImsonFZw+DVM14xOUNIsomKIYhgeHlZHjhzp+HF7GacDz1Sl0TYBROUUDVNjpp3YIjr6xGTDjSObEWTgfUPJ5yyMrVuGuzXZoZct3aq1sysFVF7ZiQ2rCg2dp6qFv7wbYRCiQ0SOKqWG/fajKaYHcNvUdaLuZQZohkI+lwhRt89rZGUBl81vfECtzCjfp4TBeX21csZhsH0RG1YVPDNHCYkammJ6AFMYY1YEM0oZzQDdGM4omCs+5q5P3mw2rG260mX1+mG3FaSYk05CYe8BzGGMSuv4tGOh1+w41JS4O5NrOo1CtSa8lZVaIpAd9SMCX2enDnulrvNlvBfA2xw0Bp2QqGjJFCMiXxSREyIyIyK+dh8SDyYTQkbEM1W9GbNMEuzq5RnVkN1ZKleaEnV3NIodv/7ajlvx3NabA4/DGHTSSVq1sf8cwHoAP4lgLqRFTJ2TRtcugZVtXFra9bpN4u5ugq2LBbcyUqtLs2DA6piot1Ix0U0+Z9XOMZ+zGursRLHSZgw66SQtmWKUUi8DgIQo90rag1fnJADGZbSfmcCdoq4rRWCnzdcdr800cwNZMGDhQnmmIUJlbN2ytppJGINOOk3HbOwicheAuwBgaGioU4ftekxC6sarXgzgHc4X1EzgNZek9zfNWVlsu20ZgMb4cABY+cAzNeeqHd7YitjL7A9j0Ekc+Aq7iPwIwBWal+5TSv110AMppR4F8ChQjWMPPMMexmsV7haKVurFZERw7dYDniLkN5ek2pB14up+AnFXXJwqlTH6hH+/0PGJork3nsA3I5eQduFrY1dKfUop9VHNT2BRJ83htwp34lUvxs++W1GqrmaMzubuN5egdWY6SSGfw2s7bq1leLp9D0D1vHRAhC3LAAANL0lEQVRldMszClv2TmrfA8zd6LhCIUmE4Y4JJswq3NQ5yTY16OKvneV4bUw2d9NcilMlXLv1gLHKYVzY56570rh7zzFs3nMMBZ/Sw3Yil33Ds8snOLN3L+vI2RASjlbDHT8vIm8C+OcADojIwWimRYBwVRvdESzOiA7da49sXGEM/9OJuNeqX6HxBmEjAO5c3X6fyoCVQT7XGM2ie9Kwpxom1r5UrmD34dOe2buEJIVWo2K+B+B7Ec2FuPBbhbvxarKge82UXXplPtfgKL1p6SLsPnw6tOnhynyu1pziscOnQ77bmyBOziht/5Ry0i3QFNNBgka42LSj1ZlzDpfnrLoMTWCuRZ7bfLHvaDG0sDlvQsPXLMRfvXDauLJvBrvpB2B2cvZipydCKOwdIkyEi5MoW53peqPaCUZT02XPFnl2c4igJoiCJhyyDRV5fePwm6nvoot3T0JGLSFBobB3iKB9SYMQduVvM7b/RMMcyjMKA/19mLj/lto2r0bVOSvrK5LOioo27QyH9BrbeWPRNfF2Y4p3txt9u81iSXMaEwJQ2DtGVH1JTSt/Z8SGTuzHJ4qYKumrG7rnYFdHdGNXS7RX7vmchfculH0ja8YnisYa8DoWDFgY6O+rmYtEqk1ATGP4hXN6NfFe/I9zOPzqOVSUQlYEG1YVtPHuQNWc5L6h/vFLgU6JkI5CYe8QJltv2BoippW/07GpM/PoYt91cxifKBrb09miaq/c7eYTOuybhX0jCirq9oo5SJKUvX+YdH23yDvnVlEKuw+fxmOHTzeYktzvtaGwkyTCRhsdIqo+lqYVvls23YlMXk8GzjnsPHgqUHs6e3xTYpICsGbHIWx/utH8Y7NgwMKdq4cCN6HwCulsBr9QSK8CaYQkGa7YO0RUES5hojycYm5632B/tmlb+FtTJeQHzBmnXvMUoM6uD8yZSe7ec8z4+UTpTPY7V9ZRJ90Khb2DRCFKuigPk0PQaWIZXbsEW56YRMW1Gr94aQbjE8XavMLcOBSa70rkNkE1GzXUrCPZnoPfuSa1Bg4hXtAU02XozBF3rB5qMPMIgJuWLqrb5hZ1oBoV4zTZRNXz1AudCSpMXRwbdy/XsOaTIOfKOuqkG+GKvQvR1Ujfd/TNun0UgH1Hixi+ZiEAYMveSeN4b02V6la++QEL8/oyxiiaZnDHykdRnbLVENKRlQUcef2sMaOWddRJt0Jh73LmTBgzDa+VyhWM7T+BDy7NeEalXJ6z6swg56bLyFlZ5HNWZOLujpV300zUUBQhpM+ePKMV9awIG1CTroWmmDZialUXJX4NLqZKZc/XBYAItCtfkWrruyAU8jn82qP+uJ/YNhM1FKZIWth5zShFUSddC4W9TbRq/w1KK849AXDH6iFMGRygU9NlXDbf/6HOKcCFJsV2ZGUBG1YVkJ1N5XQnC+mIIoQ0ipsDIUmDwt4mmnEGNoOXAOWsLBYYwhGzInh44wo8OLLcU9xMou9kXt/cf6NmxbbqJyjWJQvtO1r0vBFGEdceVX4BIUmCNvY2EVUJAT9MRa4WDFi1mie6bE2nAHqVBzaV9nWGWDqrLALAfCtTGyto/1BdIlMQR2irIaTtqKBJSNxQ2NtEVCUEvLAjWZz1W3Sp8IC3cPmJW5C4eaej1rnvB5canbq68zDFw3cijjzKpCdCkgCFvU2EbZIRFndCj12/pdlsTdM+OtE3JfXoImiCrLqD1rEhhASDwt4m2v2I30wMd7NZmm7RX7PjUKjmFX6r7qB1bAghwaCwt5F2PuKHseGPTxSx/ekTdeYOd1PnMDcd09PIfCujNalkROrKFrgxPQXkcxZNJIQ0AaNiupSgYXq2yUYnuM1WMjRFo2y7bZk2Rb+ilOf4o2uXNMTLWxnB2LplgeZDCKmHK/YuJagN3y+BySZsJUOvp5EteycbMl19x3fnQbEzESFNwxV7wgiarRokhnt8ohipLTwIIysLmDGULzCNv/PgqbqG2gBQrqjIY/4J6RW4Yk8QYUvXeq2a7bHCEFUESthQz07F/BPSK3DFniCizFYNaoKxiTIUM2w2J9P6CYkWrtgTRNCVa5CwRa/V7iMbVwBoXyhm2FDPdsf8E9JrUNgTRBATRlBzjWmsQj5X2y8KITfdZMKEenZzWr9AoDSFf4XeXxIjFPYEEWTlGjQxqROr4Gbb2eno1rT+PulDWTWGkvYJv1okPmhjTxBBIl2CmmuiqHzoR6cqWCYZnah7bSekE3BZkTD8Vq5hIk7avQpmNAshyYQr9gTiFcuepPrhjGYB8vPyobYT0gko7AnDr/NS1CaWVtr3JekmExdrF68NtZ2QTiDKo8lxuxgeHlZHjhzp+HG7AVPlxEI+h+e23hzpsdzOT6CxCUeQMboxmiUqfve7v4upD6Yatufn5fHTL/00hhmRNCMiR5VSw377tWRjF5GdAG4DcBHAKwD+nVKq8X85CUy77NY6AW6m9K+bbo1miQqdqHttJ6QTtGqK+SGAjyqlrgfwSwD3tD6l3qYddmuTecdUR4bOT0K6m5aEXSn1jFLq0uyfhwFc1fqUept22K1NK/Os6JNoesn5SUgaidJ5+gcA/sb0oojcJSJHROTImTNnIjxsenD3MAWiiT83rcDtdnpOes35SUga8RV2EfmRiPxc8/M5xz73AbgEYLdpHKXUo0qpYaXU8KJFi6KZfYpwmksA7x6mYTGtwO2bRjuTmNIOwx1JEvF1niqlPuX1uoh8BcBnAXxSxRFikxKicGSa8Cov0OvOz1Yx/ZfnV4HESatRMZ8G8DUAv6eUmo5mSr1JO7M4u7nIVtJ57+J7obYT0glaLSnwZwDmAfihVG3Ch5VSX215Vj1I2OYUYeHKvD1cMXgF3j7/tnY7IXHRalTMP1FKXa2UWjH7Q1FvkqRmcbaSmdoLbLpxE6yMVbfNyljYdOOmmGZECIuAJYYkmkuiLMubZtz2dNrXSdywpAAx0snyBt3KLU/eojXFfHjww3jmC8/EMCOSZoKWFGARMGKEZXn9eef8O6G2E9IJKOzECMvy+mNyktJ5SuKEwk6MJNWhmyQ23bgJ87Pz67bNz86n85TECp2nxEgSHbpJ49bfvhUAsOvFXXjn/Du4YvAKbLpxU207IXFA5ykhhHQJdJ4SQkiPQmEnhJCUQWEnhJCUQWEnhJCUQWEnhJCU0ZXhjrrGzAzBI3Fx4NUDDHckiaLrhJ2FqUiSOPDqAYw9P4YLlQsAgLfPv42x58cAgOJOYqPrTDFenYYI6TS7XtxVE3WbC5UL2PXirphmREgXCjsLU5EkwSJgJIl0nbCzMBVJEiwCRpJI1wk7C1ORJMEiYCSJdJ3zlIWpSJJgETCSRFgEjBBCugQWASOEkB6Fwk4IISmDwk4IISmDwk4IISmDwk4IISmDwk4IISmDwk4IISmDwk4IISkjlgQlETkD4HXNSx8C8Pcdnk5c9Mq59sp5AjzXtJKkc71GKbXIb6dYhN2EiBwJklWVBnrlXHvlPAGea1rpxnOlKYYQQlIGhZ0QQlJG0oT90bgn0EF65Vx75TwBnmta6bpzTZSNnRBCSOskbcVOCCGkRRIn7CLyJyLykogcE5FnROTKuOfUDkRkp4icnD3X74lIPu45tQsR+aKInBCRGRHpquiCoIjIp0XklIj8SkS2xj2fdiEi3xaRvxORn8c9l3YiIleLyLMi8vLs/92uaomVOGEHsFMpdb1SagWA7wO4P+4JtYkfAvioUup6AL8EcE/M82knPwewHsBP4p5IOxCRLIA/B/AZAL8D4Msi8jvxzqptfAfAp+OeRAe4BGCLUuqfAlgN4I+66ZomTtiVUu85/hwEkEongFLqGaXUpdk/DwO4Ks75tBOl1MtKqVNxz6ONfAzAr5RSryqlLgL4LoDPxTyntqCU+gmAs3HPo90opd5WSr04+/s/AHgZQNf030xkz1MR+QaAfwvgXQA3xTydTvAHAPbEPQnSNAUAbzj+fhPAx2OaC4kYEVkMYCWAF+KdSXBiEXYR+RGAKzQv3aeU+mul1H0A7hORewD8BwDbOjrBiPA7z9l97kP1sW93J+cWNUHONcWIZlsqnzR7DRG5DMA+AJtd1oREE4uwK6U+FXDXvwJwAF0q7H7nKSJfAfBZAJ9UXR53GuKappE3AVzt+PsqAG/FNBcSESJioSrqu5VST8U9nzAkzsYuIh9x/LkOwMm45tJOROTTAL4GYJ1Sajru+ZCW+BmAj4jItSLSD+BLAPbHPCfSAiIiAL4F4GWl1J/GPZ+wJC5BSUT2AVgCYAbVCpBfVUoV451V9IjIrwDMA/Cb2U2HlVJfjXFKbUNEPg/gvwFYBGAKwDGl1Np4ZxUtIvL7AB4BkAXwbaXUN2KeUlsQkccB/EtUKx7+PwDblFLfinVSbUBE/gWAnwI4jqoWAcC9SqkfxDer4CRO2AkhhLRG4kwxhBBCWoPCTgghKYPCTgghKYPCTgghKYPCTgghKYPCTgghKYPCTgghKYPCTgghKeP/A6atC8vNDYQVAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "nPoints = len(data)\n",
    "\n",
    "# Plot the original data in blue\n",
    "plt.scatter(data[:,0], data[:,1])\n",
    "\n",
    "#Plot the projection along the first component in orange\n",
    "plt.scatter(data[:,0], np.zeros(nPoints))\n",
    "\n",
    "#Plot the projection along the second component in green\n",
    "plt.scatter(np.zeros(nPoints), data[:,1])\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## PCA as a strategy to plot complex data\n",
    "\n",
    "The next chart shows a sample diagram displaying a dataset of pictures of cats and dogs. Raw pictures are composed of hundreds or even thousands of features. However, PCA allows us to reduce that many features to only two. In that reduced space of uncorrelated variables, we can easily separate cats and dogs. \n",
    "\n",
    "<img src = 'catdog.png'>\n",
    "\n",
    "You will learn how to generate a chart like this with word vectors in this week's programming assignment."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": false,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "nbTranslate": {
   "displayLangs": [
    "*"
   ],
   "hotkey": "alt-t",
   "langInMainMenu": true,
   "sourceLang": "en",
   "targetLang": "fr",
   "useGoogleTranslate": true
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
